{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Contents" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This notebook analyses lag-frequency spectrums of the light curves simulated through impulse response approach. First, a simple case with delta impulse response is covered. Subsequently, an energy-dependent impulse response scenario is analysed." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setup" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Import some useful libraries." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from matplotlib import pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Import relevant stingray libraries." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "from stingray import Lightcurve, Crossspectrum, sampledata, AveragedCrossspectrum\n", "from stingray.simulator import simulator, models" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Initializing" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Instantiate a simulator object and define a variability signal." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "var = sampledata.sample_data()\n", "\n", "# Beware: set tstart here, or nothing will work!\n", "sim = simulator.Simulator(N=1024, mean=0.5, dt=0.125, rms=0.4, tstart=var.tstart)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For ease of analysis, define a simple delta impulse response with width 1. Here, `start` parameter refers to the lag delay, which we will soon see." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "delay = 10\n", "s_ir = sim.simple_ir(start=delay, width=1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Finally, simulate a `filtered` light curve. Here, filtered means that the initial lag delay portion is cut." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAEFCAYAAAABjYvXAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAACIxUlEQVR4nO1dd7gdRdn/zbktvVcSQgqhhBYghF6kJnQQFRsIKKKgIqAG/BREUVARQZrSQaUpAhI6BAhSA4SQ3knvyb1pt5xz5vtjZnbfmZ2Z3XPuubltf89znnPO7OzM7O7svP0dxjlHihQpUqRo38g09wBSpEiRIkXzIyUGKVKkSJEiJQYpUqRIkSIlBilSpEiRAikxSJEiRYoUAMqbewDFok+fPnzo0KHNPYwUKVKkaDX46KOP1nHO+9qOtVpiMHToUEyZMqW5h5EiRYoUrQaMsc9dx1I1UYoUKVKkSIlBihQpUqRIiUGKFClSpEBKDFKkSJEiBVJikCJFihQpkBKDFClSpEiBlBikSJEiRQqkxKDNIpvL44kPlyKXT1OUp0iRIh4pMWijePCdxfjpv6fhsQ+XNPdQUqRI0QqQEoM2ivVb6wEAm7Y1NPNIUqRI0RqQEoM2inQDuxQpUhSClBi0cTDW3CNIkSJFa0BKDFKkSJEiRUoM2hK21WeD3xypnihFihTJkRKDNoJnpi7HqF++hFkra7RyhlRPlCJFinikxKCNYNLsNQAQIQYpUqRIkQQpMWhjSA3GKVKkKAYpMWgjUBaCVC3UfvHM1OXYUpeNr5gihQUpMWgjUHEFqWTQPvHZsmr86LGpuOapz5p7KClaKVJi0EbQ3nyHrnt2Bva59qXmHkaLwVbpSbaqpraZR5KitaK8uQeQIkUxePCdxc09hBQp2hRSyaCNgJv5J1K1UQoDm7bV4825a5t7GClaKFJi0MbA5OofGpRbD5Zu2IahEyZixorq5h5Km8R3Hp6C8+//ADW1afLCFFGkxKCNoC3YDF6ZuRoA8OSUZc08kpaHTdvqkc3lG9XG/DVbAADZXFuYLSlKjZQYtBUotVDzjiJFE6Aum8Po61/BL5+d4ayTZqlN0VikxKCNoTXbCHbUevbfT1dg+vLWo4qqywqJ4L9TV5SkvYh9KUUKpMSgzSBNTJccP3j0E5z6l7ebexgFgz7hT5ZsxNAJEzF39WYAyZgA056UIgVFSgxSpGjhsK3zz01bCQB4c05y7yDVTioYpLAhJQZtBOoFzyjur5Fv/N1vLsBbzeSG2JpVXS0ZpbqvDbk8TrrlLbwxZ01pGkzRIpASgzYCXmID8o0vzMZ5939QotbaPhas3YLXZ69u7mEkQmNViqtrajFn9Wb8/D/TSzSiFC0BaQRyG4F6wRX3x1L2eofiuJvfBAAsvvGUkrdtW7oblYuqkWqiNA9W20QqGbQR5IMXPH1D2xpKp+MXcyPfyPbyckCZlBq0KSQmBoyxMsbYJ4yx5+T/YYyx9xlj8xljjzPGKmV5lfw/Xx4fStq4WpbPYYydRMrHybL5jLEJJby+dgOTW9vR7oMNjQyIAuxjnrWyBkMnTMTSDdsa3X6rhbwt9P4Uo+oJ5kYjRYN8YJ9qVDMpWhgKkQx+BGAW+X8TgFs457sC2AjgIll+EYCNsvwWWQ+MsVEAzgWwF4BxAO6UBKYMwB0AxgMYBeCrsm6KgiDVRM3Q85TFGzDy5y/gnfnrSt724x8uBRBGJzcnlm3chrpsbof3G6oAo0+3EHUg5/p3sShGMthSl0X1tjQNRktGImLAGBsM4BQA98r/DMCxAP4lqzwE4Ez5+wz5H/L4cbL+GQAe45zXcc4XAZgPYKz8zOecL+Sc1wN4TNZNUQBCyYBp3zsC7y1cDwD434LSE4OWgtqGHI64aRJ+8uS0Hd43t0kGjVjQGyszqnEUMsUO/e1r2O/6lxvZc4qmRFLJ4M8AfgpA6QJ6A9jEOVfbKi0DMEj+HgRgKQDI49WyflBunOMqj4AxdjFjbApjbMratWn2RYrmTEwXBDOVSDPVEndrU1HAaq/pYjB/zRYs37S94PNKpfBTi3e+kUYDpREsRDLYnO7A1uIRSwwYY6cCWMM5/2gHjMcLzvnfOOdjOOdj+vbt29zDaVFQ3FqmGV0CmtJK0WLipBpBp47/05s4/MbXCz4vif1nR5qIlJqoLDUatCkkcS09HMDpjLGTAXQA0A3ArQB6MMbKJfc/GMByWX85gJ0BLGOMlQPoDmA9KVeg57jKUyREU+6BzDnHonVbMbxvl5K3rffTpM03Ds04Nm58UxTlWVoim0Hqvty2EMtHcs6v5pwP5pwPhTAAv845/zqASQDOkdXOB/CM/P2s/A95/HUuWJtnAZwrvY2GARgJ4AMAHwIYKb2TKmUfz5bk6toRuENPVIo17Ompy3HszW/GbozS2EXmttfmAQCWbNjauIbaGJLcV9e6vK0+i231uoqmsd5EYbR7o5pJ0cLQGKXCzwBcwRibD2ETuE+W3wegtyy/AsAEAOCczwDwBICZAF4EcCnnPCcli8sAvAThrfSErNvs4Jzjo883NvcwEsGkBaV8Tz9bVgMAmCeTopkoFYOo9Mqra+pK3nZjoRbQ5vCtty3eSY24o375Eva9TjfcNjbOIJdv+XEGE6etTDfxKRAFEQPO+Ruc81Pl74Wc87Gc810551/inNfJ8lr5f1d5fCE5/wbO+QjO+e6c8xdI+fOc893ksRtKdXGNxd/f+xxfvOsdvDar+d0a48ATiO5PTlmKo/8waUcNqWi0xAyszRp163EJTTKcrLH6NzYGJXAtbaGiwfw1W3DpPz/GT578tLmH0qqQpqPwQO0M1RoCntT7nc3lnf7+P/lX07pFtsRFvFRoTm8tazoK838Bt76xT6mlB51trxexIMV4brVnpMSgjUAtxDe+MBvz1mzBEbv2EeWWN59zXlLjX2C0LhEtoGNuKUblJJJXUyFfhF+/D429p4HnWgtWE6UoHGluogRoIeuRF+oFnyelmY3b6mPrJm475g60hzWhWSUDz+0PgwwLarFR4wltBo1qJkULQ0oM2gjMBSPjWSTiloKmymt03bMzMOY3r8bWo92XmtAUe23NaTMIXEtL9FgabUCWA/lwcetwrkiRDKmaKAFaAwNkcu8+rk0siO4KTaWaefCdxQWfU+qxFNseb0bZwEbAGpWOorFxBo3PSdikaMu2q6ZEKhkkQGuYWuYL7tNtx11PvsDVIthOsaCzCkOpluBix1ioZDB16Sbc9OLsInuz921DmKW2gPYKvAu1DTmM+uWLeOEzsdVmtsTU4I8vzcFTHy8raZtAy0xr0pKREoM2AtfrbTcgF9ZW3EvVXGmzi0GhhE6h0NPOvON/uOuNBUX15RwDeTKN4X4LvZaV1bXYVp8LiFs2V9rnfPuk+bjiidQNtLmREoO2AuP9nLp0k6eq/2UudsFsChSTIdPfXnyd2oZcZH+GMG1zacbhw5tz12LB2i3Bfx9BD4IMCxDPiiaI8tuMW2ipSNVFhSElBh60ptwrhUx8uhZUb2uIcPTNRQt6da4EAOy3c/dwLPI7yZOo3tYQGyCYZCHc4xcv4uw737GetyNUD+ff/0GwjSYQ82yVo0AB7Rf6fAM6I8/LtXBikKqHikNKDBKgBTHKTphj7Nu1KvacZRu3Yb/rX8Y9kxdq5bSt7z4yxbkYfbJkI7KEg27sfTp5nwEAgD0Hdou0mYQw/+CxT3DRQ1Owuqa2cQMB8Nnyau1/sd5EpVCd+Talacjmkc/zIBq4UK7/llfmYr9fJdtnQF17rjW8ECXEhq31jU773RrQbonB5+u3tgodNwA8N20Fzrjjf9465iLg9yYS30s3iAjNV2fpOfppWy/NCDltuiB/unQTzrrzHdz62rwdwoklWYQXrROqlboGt4GzsTaDQq+0FFPM1sTmWpHH6frnZuLHT0wNxpXk+miVW1+bh+rtheXwyZZgi9PWgnVb6nDAr1/BLa/Obe6hNDnaJTGYtbIGR//hDQy7+nlvvWAzkGYmGpf98xN86rEBANEFQzEy1iRnsizYTjFybrzaSHHfs1bWOMdQLP43fx2GTpiIdVvqClJ/KScXH+ForGtpoarDUtwTWwTys5+uCH4/M3VFMC4fA1sIwbAhSHnSBrjk9VvqsJDYZVxYt0UkTXx5RsvPT9ZYtEtisLI6Wc6SlqZ79Eky5iHf+x51Q9X/u951Wo0uiqU2rShp5KPPNxbNkbtQ7DJW7PpXSjWRD0oS9M4R4zsp6POdtmwT7nlrobtyM+HRD5bgt8+LLdqTMBDH/OENHEvsMoXikyUbC362P3nyUwydMLHoPpsa7ZIYVJaVJaoXvmBNOBgLXpy+EnMt6aK9C3ykrtvgySM/YhqzYLNMD/z5+jCJX+Nz3uj/c3nyWiegOEk8j4rniovzairN1JGSnNeOzBL3VyyBWl1Ti9Nv/1+Q8kRh9qoaDJ0w0evB1tS4+qnP8DeDSPmYuaTbcNpu1aQ5a3DWne/gkfc+L2iMT360TLbZMiWrdkkMkqK51ESX/P1jnHjLW5Fy7yiMMfrGrCZj6KmjvzSuc+lC+ML0VQAQWRgUVlXXFrw4+K4vyRocXI8v4K7IR6kkg8INyMX1V2wbiWwGRY5D7QNt4o05YtOj52VQWluHymI8b3W8mqk1oV0Sg4aEEZSFcFvF4sIHP8S/P0oWfZlEBaDgtxmo9sR3VE1k2gyibUw3vG1MHPPHSTgzxugdGZepvrKUJTk/PhWH/v9vby3A2s11jjOCmnJMhdoMSqAmSlJHXleSqV0oZxpXvVzecDM2o7lB7/32+hxmr6rx1N5xaKGCQfskBt06VCSqZ/pXNwVen70GV5ZgE46ozSC5fcEkBi7CQqvF3ZNaj0dPYRAdbUkg1ieJBTDHPXNlDX77/Gz86LFPYtoW380pGSRpyicZFDuf41xJFTFoyfEHP3j0E4z78+TIFqDFoNHq0EaPoGnQLolBUrVPKBm0jMfnG4U5Ru8leqQGIHp/bPeLnlu6+6Sfn8kwLFwr9kO+8YX4PD+hmshdx7wWtYApV01n2yU2ZBeCJPfVlPa8dQrs/5F3/brxsjKxjLQULyOVsZfeiw8WrQcA1DtUXXEopZPEyurtuPGF2S0udqFdEoP1W+JUAgLFJAFrShTiIeS1GRgGSZOTjrYl65E3gs7jppKgOlSUJZIIFJIs2MV7EykDcnNkLU1et9A4gySIVaHtwFQdhaCp1trGToHLH5uKu99cgKnLNpVkPKVCuyQGl/z940T1Wtjc9nKIrgU8Sd2kNgNab+zQXgCAYX06uzsqEFFvonyBqgfdMG4DvbZZK2uwuqZOnuPvJ6lkcNWTn2quyyU1ICd4psmIQYGDirnovIOpaC7YEicyi7SQBL68UMVCSaHlLYx6tkti4MKWuiyGTpgYeEUE3kQtTJyzIarn17l/6zFXWw4pg07dg4cLYnD4rr0LHGlyZHO8IGLgS9sQVgp/jr91Mr7z8JREbecTUoN/fbQMh/7uddJd8XNHea3Y2th/SA+jJJ4QhkFnhY0jbskqdTLBxqKxROnZT1c49xHX+imymznSbbyirGUtvy1rNM2MxeuEfvr21+cDCCdVSyEFfjWRwc372nGcoxC1GbjHUkpu0BxONs9xzO59AQDj9hoQf36CPlwLYVJur6nXu9qGXPD7SpnW2TY2ylXSRSnJdRRKoOJUY6q1pt4T+d7JCzHPEn9Tavzw0U/wtXvfj63XWAmhLEYymLmiBhu2urevLTVSYmDBzJU1OPbmNxIFneXyvEUEkRTjTRQaXP02A7uoHHKDTbWfQTbP0b2j8Pwa0S9eHRVKPPH2kkJRrM2g0FtC8wRtrsti7ea6sA3Hop9hLPi/SDI01rFEfiRDnDZjRwjOnHP8ZuIsnGCJv3GfE/4O5miJxkLbLBZxrrgn3zYZ59z1jrdOKdGuicHOvTo6jy1cuzV42nQBWbphG7bXC+5tdU0tRlzzPB77cKl27vw1m3HtM9N3qHrJXOR86hWzrjmnnYSFzH66PpUuVYTecTaXT5RzJzg/gZqoWHpVrDdRY2bArJU1OOiGV8P7Qhqj0luhuudCx5RUTdSUkkFB6kIj91ZTobFX69skSN3ThR7iXmq0a2IQu4OX/KYLyJG/n4TzH/gAAAK3x/98slw778IHp+Chdz/Hkg3bUEoU5k2kym1svfFtIGJA9vTXlN412RwvyKMriWTiMrDGnWpLFpcEpQjw8j1CQOie1f/+3dypyxubqM6FQjb+eeLDpUWlGC8kbbbyfqIqFu57HzywJSh0SdSFwrd9qCvauynRLolBny5iE5W4lyL0Vxb1lJvjB4s2iHJH1k+FlmBQs+r75bd6WbYa7puuOAMt6MzSV2OXGJvNIDwW33qS4KyiJYMC6laWl/a1svVNn2tZJlQTJQkCLrVW0yI4OvHTf0/Dtx74sOA+CpEM7v/fYgDAui2UGBR30UmJ8xf++Ab2vz7ZvhAKipm0IS62oynQLonBLr2F/jlufph6xjmOcPYdtej7deGOcs9kvn2SMJRP+Xyjdtx873zv4YeLN3j11I1BNp8vKPo2zjAujvnPdZ+XXBVCA5sKXYKs6UMsthB6jRVl4ZhyHm4zRiCMPS/ueFI10ZoiJAM6B+NSodhQ7LXbAy4F6OUuWrcVG7cVti/ET/41zXnsBpmBdUeifIf32IIQxy2YaiLnQmKUl1o95OqHwqn+sOYmkouLs0GHZKAZMEXZjBU1mLGiNDlfzNEo2wyQzL03kfRQtAFZfDe1AdnahnU81GaQCa7dFQV895sLyJgKG9RRI/vimakrnMfVWN6TEnMcilGv0Os99S9vAwAW/PbkWI+cAEU8hxHXPI8TR/WPNhVIyi1A9C8h2qVkoLCiuhYPvbPYedzMWhp5z+R/F0dU7vAjvnfyQsx3ZPukL+rJt07WuO6aWg/nkZBQAcBTHy9HPu9eFl3XGdeuz+soCcyqf31rYXB+Ep3xVkk8fKkl4miK090271cJuhsssHrC+0oFgHIiGbiIJk3nUSiB6lzlT/mu2ovbgCmsX/jKbDtlxDXPawyDgjcCvYCuc3keZOfd0bjw8GEAgAMi8SRNh3ZNDADg2mdnBNv4mRPFzFpqTmJX8rLecmN3G9OSywsXudNvf9s6Hvouz1xZg1teCbfbu8mTnycZjy/wh5fm4Ompy9GjU6X1HFM/m3QzeLuKw3uK9/wzRw8ihvDk7Yy/dbK7D6cBWUlL9vMUMco0wxtjCyCkwyzPsOC/Lz+Qeq6FSkfFElAX1m+tx4vTC0x37ehis49Bij+9KBRiIykGddkcunQQSpuhJYzwj0O7JAbOBcExZVR1+lLMWllDPA3MduL73mbhaEQf+tmriH51e4P9HNpu0vKN2xpwzgGDAABfGbOzdwxqEdHURAlfr8a8hF07lBPDaGle59iFzVUuD0xfXoNH3l3sbYOm6Ch04bXVbrB4lug2g/A19jlFNDiYHkBXI0X7ch6SffqP21DoNpJuNWgUx+/ZDwDQqTKUaGx2l0KgO08UKSUmxJE3TcJtr80DAJTtQC+UdkkMTKjpYU5qcwLRiTT+1snOKFzlpRPnEmhLp+t78fyxA45yp/oovJqKcmNzG2PtsaWwbpIQCvP+I7znLjVR4SoHf33XokPLf/HMDG8bPTuFKdILHZ5NzVNvcRHSbAZlmaAfn2TgIwa+rLBJ8zYVAtNpIQ6FdDGwu4gfGjWwW+z5xcQCNbVksIYkBmzqqG6KlBiAGoi5o1z/Do7Lb9fzsk0z+hKP+uVL3uNmI15i4Fr0Pa+RehFM7iO68Ba/8jc2Klldsuul/W2BXhdx776bGCTvgzt+JzrXcoIt7TIdT0UZC3z3fbe7EM+suDHpYyn8GRe6EU4h88iWeys43WhmYhG7s4XvfdMv1NNXVKO6QC+lYhFLDBhjHRhjHzDGPmWMzWCM/UqWD2OMvc8Ym88Ye5wxVinLq+T/+fL4UNLW1bJ8DmPsJFI+TpbNZ4xNaILr9MJlIFaLYrCLlMv/voBJUeh7Q3XU3kykMSouE9k8hwqANMcftRmoev524yShOETqch402uCI1rxn8qICevBJSv7jhSx4mm5f/kmaR9/Wj+3azQjkK56YGj8uTx+FjsnWbiEY3rcwXXghkq+truv9SBrcRc8O1KaJznTjuD36xdaZsaIG597zXiN7SoYkkkEdgGM55/sBGA1gHGPsEAA3AbiFc74rgI0ALpL1LwKwUZbfIuuBMTYKwLkA9gIwDsCdjLEyxlgZgDsAjAcwCsBXZd1mgLEIKoOb4k4dDLNrUhTDFZsvXleyK5t/X2NHuaP+mpq6YHyme16hm9v4+irkFtgin9U99/nPFyLqF6vyKFbCUalKkuqqbfe6PpeTbdDxhL/LM5mE6Tr0+VwqFHNvjhzZt8A+CqnreU+M/zRGw3curXXzy3MAAIvXNy6+pmdnuwOHiVkrd8x2nbHEgAsoP8gK+eEAjgXwL1n+EIAz5e8z5H/I48cxwXqeAeAxznkd53wRgPkAxsrPfM75Qs55PYDHZN0mg/m4XemeA8lA1cubi5WSDPTz1ARLyj3rY9H/n0D8nIshBq5zVmzaHnA4pteTeZ22fPVJX845q5JnmYxoyHh4j3268JdnRt3/ah3GdhdNUf041UQFaDVoC+r6k94v22U2ZHUJVdQrfAEOT7Gf61LdxHVlk4RKjca4KAPubVijjFB8+50rhadPhwq/y20cmuhWFY1ENgPJwU8FsAbAKwAWANjEOVcW0GUABsnfgwAsBQB5vBpAb1punOMq32Fwcf7q3Qh1rQYxkH+31JrpHDx9xXKm5opI2i0mXYmju7HDehF3Sb/NoDGeGKfd/nYk4nTohImxQy3LMHDw0DDqSeplE/U/chgoY++/o7ygxZfULfSO2RY9ZUCm8yqv9VGY1OGan04C6rj2XXp3ihxPuvVlqRL+mdc+dMLE4H7Z7qWpAjWdP5J4rZ07Vnjf7UkM1MWgpWynq5CIGHDOc5zz0QAGQ3DyezTloFxgjF3MGJvCGJuydu3akrWrHkk0jkBfBDtWlBvHxXc0nUPhHLzZZji2sMAXeEX7rCSuhq4zBvfs6DQgO7nnGA7QNbxlm7bbDxig18Bke+p+FLq/bt+u9oRtcTaDJN5EcdDVOX6JI9pPtMzGsdNnVKjU4arvcl121bdtluMj2o1BIbQ4MKZbjh124+vBviVA9J0yn5Mv8V9jNyrL5XkQ49QSUJA3Eed8E4BJAA4F0IMxplbHwQBU6s7lAHYGAHm8O4D1tNw4x1Vu6/9vnPMxnPMxffsWpnP0wWlANmwGSfWLPmNknBeFy6MJ8OvGVb3K8oyR995+Toax4HpNUdn1guiLnHUU1r66VCXMekJOZwy6a6lHLLKlJIhb9F1dc9JNXZakwyhyjVOn0fOnLHanbbARjZmWdB/0uRYsfbiMqQ5Viqv9MCiTSgbJFreCs786RmELhFQGd9ez/pTsPexiABVsr6t6lr99fjaue9bvZuzDM1NX4MQ/J9+foamRxJuoL2Osh/zdEcAJAGZBEIVzZLXzATwjfz8r/0Mef52LO/4sgHOlt9EwACMBfADgQwAjpXdSJYSR+dkSXFtiuNRAamIoomCKkMWsD+8uWJ9oLOEYouOxnkcCYZLGA4TBZDEGZENdJvqLwuVVldRXWpcMZCbOBJLB4J6dImWuxcOdwlo+Y3L8lZlhYFSx3kRqGHTROefud53n2vqx5bSnt6NQPb3rVq537KoVZ1Slh5NKBrZF/HceN2E3EY9XBZn4i9zJ0FY3+t8f4/GgJ51NEvgyl+5oJJEMBgKYxBibBrFwv8I5fw7AzwBcwRibD2ETuE/Wvw9Ab1l+BYAJAMA5nwHgCQAzAbwI4FKpfsoCuAzASxBE5glZt8kQmVgOTl5x4i7JIe4ltE3ULXXuvDm0L1sbSTa5pzuP0fLo2MLxR9VELsmAcKOWdl1ic9LFSqvGRH+q7yQbgTjbouUxY6CLwWX//CT4XVicQbRy0vNt41a6eb29RkgGjptjSr6+MQEIqAFtj+7U5h2DZdR/fWuhs34hxFgxDq4zqA0r4kZtrP02wlIKw29iaXkHInZEnPNpAPa3lC+EsB+Y5bUAvuRo6wYAN1jKnwfwfILxNgnCLJ56eS4gAuq4W4VjbddyvFCbQVIOUB0xOS4fh1wjDd+T563Fj44fGRyLSEC267C068pdlPTdMW0G4GHfxURfx/VhayNOcigU8Rliddj679e1g6UeHZtxLM8jTgFJkIT7vurE3fDHl0W+LJtk8P6i9Yny6RSqKneNzTYtcjZRlqBP16pg7rveebN9rpU1nhq4CG9zIo1Ahtt4qCZsoCYixzOsODVRXPbNyCJL6m/Y5t4cW1XLJJUMeCjixhnAw+A7z8DhlgySvjyaYBBkjBXfDb48/QUQXVpu2/a0JAZkUjVmXbL0k6xvk7jQ/FJxxvZCd3uj5VXloTulUv/R7kb07eLtO24MhcJGZOMkg28fMTz4bb6PcWojwD/2pDEvnSpbnmSQEgOEk+ZzI4gkVBNB+wZEoE+cbtJ2NPYckxaQ30s3+LxyRE3GdP486Z4FFCbXltSA7LIZcJ6MM45IXgiJ40aHPtvZlscZUWHfQT0ixa7nQ2lR1w7+F9l2qcm9iSyLj01VYfzeqUdI2AqdY76+RXn4mz5iM8U7ANwqE6zFQY3R5c5qohCHAKVSdF1Px8pw2Xtumr5PgytJo17HPc5/fbTMfZCgqsS74ZUCLW9EzQA1AcwEZCrCMCAKZBZkMvEvnU1/GrckRI23CblqWY3Bvl9rtB/g6wcPSTYG9ZfqqS0vmstOnOc8YYQsaQsMnIdxBovX+zYMsnHOrrGEv22qFJcAktOu3TMU62gaZzOIlRY418YX59FTqD0lzpOHzoVC024kzVFEx3DuQTt7asa/l/k8MKCbUL39b77u0GHOfdu99DE2P/23e/cyhcN37d10KU8bgXZJDMxH6Xq278udm0IDclixPJPRUgfbcOYd/4uUxc0Bn2TgPY90oHsTuVUCbkJhcOgWycCGjo6IzKSSgWYzYIVIFLFVrHWp2lYVJ7EZFLc5S1KiHq1n9WghCx43zitGYvWN0dWcuS0sAHx5THShtrUbbBaU2NMs/E1die0qQu48BgjC7urWmZdLaz9utH6wFrpHWrskBibiIgGVI4vGVbJQR97PEeBkw147+aMWowtx+LurxwNB26M3gc0gz7nz5XdtbqMxo7YxOMeWjDOmVZj8n2SHM+tYEqhCyiw71Zj9NVijf2PGQ9ro07XSOUb7+KJl9RZPKkMw0J5ZsTYD12n0emyLd+y8sBSazhlx4Npzo5KvhXjGEAPOudPd2Z2uhNRxMEtJITz+Wh45SIkBEPumKi6GvnBlGUYmtPvc9VvqtP9x08acV4VGrmYikoGjH5gLipuzDI2gflWJW/3AE12Hrptm2iLnc8UrVjIo1zhM+SyNa1dRuYWmf1BZOVUMRDE2g9P32wlAvK87h64mCoMlPQ+/gAN0PtC5xQIDcnjcpta0XXuQ6iWhVxFtgS7ktkt8c+5a2Yfretz9JGE+XMGpSWHa9VoKUmKA5Hp8kzsJ54C7hWnLqvW+YiZbsSJosEk3Y7rNwKP2yDsIgIsg0WJ71tLCOE7bmBSYbC+MJi3sxrhTeoflGY3DtI9VLVaFSAZASGjUfU16D+gzGTusFwB7jIXpWkpvj28TG6Dw5+RUE1mOWw3glvODgM7EkkH4u8zy3CiUo4UziI6cNbin7lGWZGF32tQSQtj1CjtnRyAlBnC/NHsPEiodm2tpWYYFC4tvAq0zJIM4e5kvfbRvzqlj5kRzi8r6JJ62PCRa0UR10bas7XrE8iScsVZF2gyChc17noUwJWCKbSafaJR5lBglIeiK98vHcemWc8PxiTbidjoz1URr5U5ZhauDXOWESFu8ifRsqrb+bJIBdx5zjCL4ZZPoCgG9X6baNpnNSD9WqJssYzt2B7OkaJ/EwFzsHEuN6aKmGbEYc6apoDBTJRS6UUih7wqLqIkckxtca5vmv3HaDDzj9PWV58m4J/oc1DXY3ARtLqjRtuwwnQDCNqPHxX/57ZGcon1zzFm9WTs/6XPUbBpMly60PihhNs5Tjg2uLpPYU1zlupoo2p6N07faDAokBrSaTaJznOUszQbvrn1cvv6TuJ/6YHr8tRS0T2JgQD1bJZYrmCK+qWJQx2tqs9ju2OA+buMYE7YNXpKAvq5JXEuFaiE8Sie0awxF2wy421hNQVXj6hpsKo9EEdJODi/8bePOXMbBnGcx8PXxl9fnJTonbDv8rYYXt9MZ51wjVq79OYL6ThWavb5LmlX3T09TYiEGnlmcnEiGv8tibAa2c8xOlR3GtMc47QwO6ef4Pftrz+LIkX3cA5LIMOacn8XsyVwqpMQA4WJ3xK76g1QT4/XZa7T/gLIZhP9rau05WUxVRCxX6VMTeU4OvYmSSQZ5rhOKXmTXJZdkoPVntRk4xobCJQNAXJNNTRTR6xcwFt2bKCzfJom5KxlhoVlCFQ+QNeZQHDTJRfq+2nzdXVHiADEgF2gbSJLcjzIaVptBwiCtJCpW19h0Bst9fqXD9TvPqerXHJejf40Y0ErcysT4oJwjbEhqQ2kKpMQA4YP2vWDiePhbNyAn5/wKjkAm/31nqmOJE9UZevyB3cP8N+YQrdxKAZJBPq9LBu7gNGh1qDgP370uREohv6lrqbLtuFx7TR29DxxAny66u/GnSzf5TwraDhtXnHcSA7JtLjrHmUBqco+JlFvOs+UTtDMT7mNxY4uLM9itv0iJMbK/PTUG5zyU+o25nciV2ZASdNVe7OlYWb09kUpuRyMlBiAvvDExzJeQHqc2AyC5G1u8msgcW1jgS4OgqjEjibVTXWL0pXHejsCbOMLkSwKnnesTISQyknsKJQMXZ1YgyKnllghk83mFjIK93NoFIbQHDe1pdusF5TIzhqrM2R90NZG5Q1+0vqPcxa3S7i0qmjjjus19lFvmlA9OYmCp2zkmIyhHuOj7chHZCB+gPyPOuZGzLJ4aKAO/DbY4h06VjdteMylSYkDg8h/u06VSHtdF1XcXhqHsLvHQF0SWtH5nORmGebJB0v2Y6VB8gUQubteVvCtON+xC3uir0pGXxTRU5jkPCLLXZmBtLZ7zsqajcNhLCnEn5GSMamGybcATNz5FrOKCyDi323ycah9He0m4VXoZtghku2upuyx5EkP7gms73RYXpLXFqZrIfK70/pN5qjFCXCuOm1Mmjtuzv5Mg2yST/t2iWWubAu2SGJi326YKAMIHo1IIm2qiidNWxvYVEUMLVRMBKJe6T+8CJI9ljIAWb4ZKjVt3Sznmjm+2cZpt6Af0e+7idGgdpSayeTJFuXTbWBxDoc/QZkB2GKej6jsPYeQkc6Y8dOKoAfYBmf1T6UgRA8Kau5LWmWoj25hpfdu43RKDfSG2BZ3Zdwazj1nUT0gMSLXymAjkbBwxIMd872e5Q+cTUdEVuAUpY8kIbzCOxu6vmRDtkhiYcAXAmO5vZpyB1oaT27L3FTcW+j/r2eBbwVxIbeVmfVcMg9mPMmBqdWLGQGHqVZ2pALR+maYL1/S0jliAZGOhzzB6PBJjAfui4QtC46S+OpT0fbYtvNSbyCYlRKQ8NeYY24B5OElcgiYZWI7HBZ1967Ch6NOlijBg5tj8YwbibQZxLt/ifsm65vtmOIkE5zjGaEq9SSRmZjbo6N82jqZESgxgN4QB0UllRiBfdMSwsA3Xwy1AvUD7ou02GC6uVshjIs4g3mYgchPp/YRj0OuaXK75m7ZpH5rel2tqa3XkG2PT2buC4uLKRIshbLmJosTb3mfeWBAUMkzXI7ucE1yg/SuOkNoMbC7MJmFXF+nKs+NU0TiGSK/dlkgh1mZAyob16Rw8W1t9Nzev3+Owbct4HSogBUpQI+88KaCJKH2BdVr0vrVHHUrqtcF2/S4JpdRIiQHIC+sQGW2cRhlj2kNyvezRhcs/XaKcabIAHXUsw1gBEchRbtI2hpzS28NeP66vPDcXevvkjqSj4HYCmCya2c0VKtjE76iaKJ4Qm1IPR3QhTuo+btM/03lXl7PEs3D7IuLbEMg2pmSupVq34tvDSJj97DO4e/BsbWNw27jC35kYNZFLBaRA1W5mDdN93AbT1kbHlmRuMjBnPduQberMpkBKDED0lzHqEVOfm2QSRDifODWRqY7IEzc4z6nqUPIIZONFoNdielUltRm4xmaI0nRzEdf5KoW1njNJveTJ+rWPRe/DhCugzR//EZareWFmzkyqG7cZa+mi7rrvNpVfXDStuZA6M3bS67NwGj5GAsbYyhgLPMVsY3TN1xWbwo2dtCgDblNrcu3bhBbE55EyKbNACZDuhad7EyURABlz7xVtu/5UTdSEiEjH3F4e4e7ILGAwPTjsfUVdS/1jMyewvhC4Tw4S1UHnvH3eRK5AKnNCxm2WEjc+Uwpx7pVMF2owmNlOQ27Szr1rZY4xuu7HKfsOtLYdnudeNCJqIvDEi5wJWq0sMCDTeWa5Vs6tRkzXYuiKQ3BKt3TeWyQDbXrEEYMM0wyoNhuZDVc88Wn4x6Di5r1WnL+LGNL57JMMqOSvezCR55H3R+/bwABs2mYnBrYxp8Rgh8I+MQPJwDJxOexcq4kkC5dW3yQG2WQTLZQM9InjVxPZ60UlJLV4xC1KrrEZ9omE6ghTBHdxvLbWfGMJ+wjvVWXgsWVfnHyOALqBlekLc4wx14TJRQM6J2tb3zj0Z6b6zDriEwq5j+aYbG6dcddmShaM9BUlSPY2XGm0OY8SvfqYbS8pcTWr0Puo566yn8PBY99pE74oZRsx2FFJ7VJiAPekDjxZjP+qMEmKgkJdS12bq4i23Oep06IRyO6F12UzSOJeWQgxyOfjVUxmuVowbJ4akfO5W9cf14eCeuGiQWfyGkyiTn6bi2XWotZJKFxpcyyQDPJ0Dtjvu+16XZJBnOTgqg+YkoF8HrRuzPlCMmCed84+JtpvRpN8o1JYQ6xk4L5ulzeRa4xUJaj+NwapZNDMCDklfz3TgyTJJEhqIFOIqIkoMUgogupBZ+4FwSkZJNjtydasT8USGrh990rn2k31UvCcIgsztxre4/rIaAuMq235bTRIo2rp7wwzJUZ7uy7Q8Slvloas32ZgLojBYm/LDUH6cBF9E3pUblwEsrs/IHTnVYRkwdot+hhikuIBUVWV651xvWu6Adlg1jTJwO5aas5JSugTvaOetd02T1JisANhm9S24+a+szmLOsCEzyXRBvMFpRuMJ5to+g6rrjM4otdTSD8udQUQvY+cHEsq8jJmUS+phdXSeZJMptE+ovfJZSh27XNA66g2NddFhwulC7btHWl7LolMs4nK7zhbT5RD9hMPwL6QxakATcmJET3RbybO0uq6iKa2j4LWt1sycF2/17XUYTPQHCy0OemO5HfBt89ZGmewgxHJjhmj11XH9TS23Om6pp1rLowxs8Wc2PWaZBDbnchamjAdRdZhnI4LVjLrm8cjLotUMsh40vcai44pGYSGz+jC71usfX2YY3cRlRrD+4NWM6WNnEUnHSd12tpVkkEcMRD98MhvX9CVra0k0q1VMohpw1R9UT97n5eWa2ymZBAlBvL6HPfcZzPQpRhiM9DGop/vijlxoVDJILUZ7EDEvbC2RS6iwnDMgWhqB/9Yop4RBU40g+/wevhoko17DGEd//W6FiHOw2M+Jkd74cFAI0UpbEFnSQLRon2QctjHru7Mkx8t00pdBuSyjCkZROuLccQTK8WZxnmtKe5YqTVUFds+CLQNUzJwza64rKWuSHZby5nAtdR+v5N4czGDILmTRNrLG7zeROFv145qumtpEUFnnmOFbpRTSqTEAPFqorCezn0VZzOIkQyM41QySKb6MF1L3Qt71kENknBnVnWF41ieh+eqoCwbIpIB7Pp3W5xBNEWFqw/yx+InaV5WEqLCjYXK9kIXE4FcYYmQdrnR5vI8kiQtLs6gGJuBbSnz7Z9tlpWXMZmEUPwf2F3fgzhJrIPhL+dUB7lUTj4vQI0YJ7AZ5E2GJeE76oLtUorZ2rMYpMQACdREgeSgcwCT562LtGHilZmrtBQChaqJlP6zoswdtWgiiVCZ5+b1hL99Cb5sv4Myl/oB4QtT5tnYw1xyqERBx2jjsqPSiOtZkoXb0rdvDwtXO7oaBVYim3Th1SKkLWkI7JKB+ASRqorzd7qWigpRycAvSQBmKghxQM/h42YSADFGRhiCc8fubPRlH8P5h+0S/E4qGah7Hpea3nYOoN9/nRGCVl64msjjWkrO//tFB4s+YlssDVJiALLYx9x1U09M85K7pIuXZqzGL5+ZHvy3TVw6N8xweEUMKssyiSZaxghBdo2Lc7EY9O1apdUT1+ZoPGbSq/Oion9YXxmGfc0r6YYj+uLZ+uYIX6Kd5CY9CS4h4qJoG7tzrOS3aSC13b+IHcUxPm0xsujU7IutkFJDNZF9sTfHktTo7spaSs/bb+ceztTk2j3PyDiDgHEwx2YfxE49QgnCVFU5XWjlt8vbSNmlKPRklK7r0ZlC9cwyLNlzZQCG9OpkrUefSd+uVdhvcPdEGoFSoF0SA86Bnp0qtP9AvHEuiY3ANi8/W17tPU5zj5jEQHkTVZZnEvmqZzLRoBzbeFWUbIVaQCgxSCIZ2DjU4HyjnIdeQWUZj2SgCAYQLBhJ0yyoPo8c2Tc6WG2M5KVMSDQBsdet1p+DMJqLpU+asYFemo2DNKUlVZbnoZrIJslGTiDHDxneS47VDm3zFsuKkec83G7VOi+M+8Pos9QnSxI1kamWijOUR+xw8n95hkWIq4sYOzdXIuri8oybYdPKGdChIn5Pj7IMAI9atdRol8QAMLwiErr/6ZuO2+vYJubsVZvJecklg4oMC4yAFWWZRLpDBma4TLoWI/FSVJRHI29dulaXntwssxlzVftlGZ/NAEEdtWCYbnyuawleSCnau1/K8LdVTeTgll1EwvxtLpY+acY+Pjuxou3ZpIxcnqOiTCcGrh3SQo5ZHFcENBGB0jLihmNWLqO2FmizZWq/jYBgmWOLZ0QaNG8g7lT7mNeJ4H+4eEckA/r+udREWlxBeCyTcTOJ9J3YWpe1VzL6D6K1d5Bo0G6JgY0rdHIl8ttcmBRHpbUR8+Bsh+kLZu6ZENoMMok4hIBDM2Bz97R5oADJvIlcumsg+vKZKSCcL0w+JBhqwbBFeUf7Dt18yy2Sjm2MYiy0nMuxG5wiIUB7DuyG3529j1Zu/na5Aboim33jS7bvgzinPpuPbADv1KWra5JjqohJkex0x1XqqJyYR8oDzIR2zzPKiUBx7XlnXb3czqwoQmg/R9aP2AxEn+UWxoTSlXIthXVYbkbu5yhxsQ9fW1t6dKxEj05i90TzEedNYrBjvEpFfzuuq5YFW8BRHEccCUO3GGAL9RYyYYbDK2+iyvIoF2ODKzeRXTLIB77srslua8tVJ4gqXbNVK8/nw/pisba3r7jq8w4dKm0G9oAek+NVUg4QvsCuW6XGWFmWsQbn+SSKyvKMldi4/PB97SaKg3CM3+b+2pDLB1KeOtoQs0gqoq1043EL8aiB3awEqiEn5pFNB6/GrFAmFzi1OCZ3bw1/a8GS3B1c5iLwSpIoL4sSL3Of8zhQqTSjvB4soO99h4oM7vjaAQCitgOTGfQ0WXK0X2JAfoc5b/wcsckR2CIZ4/yEfYsooE/cDAslA5t+0wbG7H7PZrd5LoKjwg08CKebJB2F7ThX39EXXA3JZWAFxMt3yj4DcfX4PYK0FVlDJQBE/ecpdxgu1u5FHQDenvAFq57IJxlkWLjYuwijafMNFiRzzI57oNsMosc5j57LubgnlQFht3PdQX1j3IF6yTG/cnmOQT064j+XHmaVqLNSReUi86a0ozLSqrYp4gI/zXNMhkE/J1pfjFcngq7+XTm+TAZFHSv3SO+6+imDvl2rcObonbz1MhmmSVFNjVhiwBjbmTE2iTE2kzE2gzH2I1neizH2CmNsnvzuKcsZY+w2xth8xtg0xtgBpK3zZf15jLHzSfmBjLHP5Dm3MZ/vVQnAOTC0d7i5vGvSmDCDf7SFSpXHGHmtnBMpM9MF12wX+kVhM/C3DYiXbeE6wplzNd4oZ5rN80DHro0hgWRgW2xVH0oEpuXUDdG3UFeWZ4KUGnnOUZ/Lo6pc51xpig6z7zKL2st2ESr4KShWBN9pM9DFdn1B0Am4DcXkJnK9Bjb7Q0MuH3jzqKNxuYnU8bh0B3kOdOtYgaryMuOeCdRn8ygvy3hVgAoqkSIlJNq1OM53SqXc4yrq6EM3IOvIOZ4lrWd6FSpBNZFzBAPOP2woKdfrmYwFlaKaGkkkgyyAKznnowAcAuBSxtgoABMAvMY5HwngNfkfAMYDGCk/FwO4CxDEA8C1AA4GMBbAtYqAyDrfIeeNa/yl+dGtYwUevOAgAISrj5mI5iRQHJGoIzmdODWRheDoE01fDNQmGOu31jk5atNtzVrHFIe5NDpmomqVZGk23GURsZ1TNY5bMshzHiy2asGotyxyNjVRVDJw9SHbh93ryqW2iEgGpA4djrl+h0xCskWPwiUZ2HbP04iB41rMQYVeMHF2Fh5IPLYxZfNcqt1c54e/TWLtivh2DBlvXHVMRPKNi4vJGS92A1ETmd3lDWbMBlMyoCpQl5SixnjdaXuhQ0WZbD/K9dO5VGZIUU2NWGLAOV/JOf9Y/t4MYBaAQQDOAPCQrPYQgDPl7zMAPMwF3gPQgzE2EMBJAF7hnG/gnG8E8AqAcfJYN875e1ysqA+TtpoUIeclF/IYsdrMSZLL57Wc56IN/4OLsylkHQv7tvocXC+KqVqiCFVcxkm8CMnA5V4XHI+OR9VVHhRdq8rd6jiujz+bFy6pVeVl8rg4L0IMiB49FP39XLHKkUP7BqKLNlUTMYdkQMdjzgc3kXGMz/Ms1XlmZk8OwZ1XGAbkOP97avikYzXx6qw1mLGiBoAuldHnUe5RNb09PwzOzKigM3VfEqrPVHn/bh0ikf/uVN2iPOJN5DMgU8Ju8ZyK/AbHXOkt6JMMArsCkcJs9gBdMyC96nYMLSjMZsAYGwpgfwDvA+jPOV8pD60CoByxBwFYSk5bJst85css5bb+L2aMTWGMTVm7dm0hQ7ciIAWKUyrAK0G5tJneOHGLvd0LhyyyDs6k3DPRTD2nb5G79rRRQXBMLs+tKgInPeP+Ok79OBQxA7p0KPdeR8iBhjEWVYZkUJeNSgYBd2YhbhTLN24P2l+2MdxK0RmVq+4fV2J7tH16Dr2fesLAZIuey/VVrxNti9oM1J1yRiAHNgXjnjn6o6htiO7BnFW2J8fC9czU5cFvIV2FvUW8iRz95oPF1DAgWw3qOpxxBhZ3bV1NRMdlZ4TyHLj5lbmyvsdtWqmSjPSrPjWRih/dQbQgOTFgjHUB8G8Al3POa+gxydE3+Zg553/jnI/hnI/p27dvo9szovdjdXM0WVme6376cQRFwZVbRsHF5Q/r09l5g+k5Px23h8G5qPHyoE3GmFRzEW8izYDs5yjFdViOB9yeEUjEObZIyaBLVbn3OtRiyhAuPOGOc6JeVDLQ3VLNsVLc+/Yi0T4D7n5zARmj+LapYNS38vum5eZ4KDGgnKeNQNrgclMN5ipPriZyu5bqx4N75qBQfbpUok8XEalOCbGqXZ/LBwbkuGtS0pUaQ1KbQeBibZG8Co1AVvPJJhm4mDGN+HuSR7ruoXpmVHiz3TEthbaSXluSZMAYq4AgBP/gnD8li1dLFQ/k9xpZvhwATTgyWJb5ygdbypsM6t6qh+HUdUfOCCEMR5xEKfKg3AerzYBy3A410e4DurknmjznF6eOwoG79LTWUS9EhdTtKj27TV8cF9EJABu31UePG32F54ULe+fKcuc9yuZCtRtjoTShjOiKYDVEJIPQo0RdT9xzMH3ygwXUsWjn8lwzOtNanyzZGPxWi6YYS9hH0kVDlU//1UkRyVD1S9sa0bdzsCBWGm61cTud0YBG85oodu7VCXsO7ArAmBuE+JdnlAHZ9q6I78N37Q0gzEgbaQ/u55bNSZuNISELRsAlAcHah+5aqp+TxICsMSPaWNxzLkxZoRN4X6K8TKaF2QykZ899AGZxzv9EDj0L4Hz5+3wAz5Dy86RX0SEAqqU66SUAJzLGekrD8YkAXpLHahhjh8i+ziNtNSlCbktOTMc9NyfM2GG95AsY9dOPky7ibLMrq2vJ+MKJU1HmFkEDDs9qcBTH1GJcVZ5BhjFs3NaAhhwP9Lm07RwXroQ3f2k/o63w90efb4QJl46Wc66l1XCt0xu3NWBzrTCYMxa2s4dciNQgbamZ1fsZZwxVqCjLGLl01BwwXk55bz9esglvz19ntRn06yryId149j7YuVeYQ6ecbOsWDcSzI88Fp96lqlwrD/X6YWqP3561D47ara8WiyLqiOPxier081yDyhOJzXZfle3JHYEsSmmglY9xsKEhn7fGkHDOPZKBuk5TMlCMQyay0CaRDGw5rQC/xx+VzIP2LfVskkFLshkcDuCbAI5ljE2Vn5MB3AjgBMbYPADHy/8A8DyAhQDmA7gHwPcBgHO+AcCvAXwoP9fLMsg698pzFgB4oQTX5gXVrat7HedFs+/g7jhqt77YpVenwGZQYehb47yJ4kLLqSrKUC+6jVOKGBjcLhWDaxvEi9+hogxgwMwV1QBC7ttMz71L707o0kFfkOLmpDpuiwNQL2BluT/h3tNTVwAQ3KNqJ9ysXtRRi9+9540JzgtE/xgud48BgrCUZRiuP2PvcIwOLtJsJ8P0+uqaAGD0kB666E8WvaQRyNRuQhcOqtfPBVymuE8mMVCIlwykTcaiKjTbUbpu5QmjxgLQCGT7dYXMSqgCdN9vt2Rgy6NF74eJsA+PAdk4Vc3N16482una26myDAO7d8Dxe/ZHngNn7S9MnLv16+ocv6mSo+O31bvg8KHoVFnuJLBNgSTeRG9zzhnnfF/O+Wj5eZ5zvp5zfhznfCTn/Hi1sEsvoks55yM45/twzqeQtu7nnO8qPw+Q8imc873lOZfxHZSMQ724gd91jO5RGelU4JTgiPRbqNJVq0UHAHbr3yX4TSfuvBvG4wfH7ir6sFyyRgwcIjht04yYpN4NddlQMnDpdml7GcZwyLDeWnncUzFf8McuPgSAILLqBax0BOaoa7vg8KEAxLVTaQIIFwpVPqxv56DffEQysA92n0HdMVBmNu1LVDouTlXNjT5dqnDcHv2CBcLmTaSicBU214Y5aHL5PAZ064BfnDpK79BAnocSoVVNxKkxVeccTfuPYgBMBF5A2ZBAq7ZtoHsljNt7ACaM3wPH79lfU/WUSTWRDx0VISGG1qSSQTYXSgYn7zOAnBDvWmraa6ibc0RNlBd2rRF9u0T2TVBYuG4rBnTvgG4dhcqzd+dKdKos83Lx9Jkp2Oqret8+crio40jx0RRotxHIAHnZ5L12Gk65WoRyQs2SidoM1PP6+dOfAQC+dvCQ4PyuHcIMqbSLirJMRPw+96DQrOJKbWDC9BdXoMY+KhlkGAvGdPFRw7X+RXti0pqcZpzuMiSa4eKoygP9dLk94Z66Lz0DVULI8ZqquAZpsAz195wY6PwLUo6HKg8tD788PzIH5N8OFRl071gR8UADvbZMlNCGBJLkXIL7XlKffl0yiObJoQZtAJHgPGogt1xSVL1krS3mejm5Z5ccPQKdq0IJQahL1aY10VYuOHwYAOCiI4fJcRPVbD6Pob074S9f3d/Ru0BDPpTCR/bvimcuPTwYW1w8hTMdhSMCWd1/X5biT5Zskju2hXPK500UBKZpaiJ3nEEZYQhajGTQFhFEA6r/8tuZ40R+hwFQQjLQ/Ktlm58s2QQAzg1tzJfF9ChwxRn4vAqU0dPc6YqOXdkMOlRkNK67fzeV/5+MMc9RxqJqBx+D0oOkBM9a1DvZfB7lGRakmTARpggIVQmBZGCofoJcOGRcecLt+cZK3WmpUEc51X5dq/DPbx8cjF1dUzklQKT9QO1QpscujB3aC2GcB9diG/xqovAeKFBHBU1NRCq5FvVT9hkY6QMIiXacZJDNR3dRU/3k8zywczDHYqiIlLKDUDVRNsfRoaIskBqSOBcAhndVwTaD8HlFXEvJ/LBlGKDIMIMpZO7xL1y7BYAIHqXX4DJgUzrVkmwGbRamz/iCtSKNg1oITKjMkMpPWqQNtr+AlJOjc/GDRRtgg6qiJvZ954/RxG7BRdjhkgyG9+kcXFxd4LNfJhZaz0KgXggfl20aOHt2qox4iFSUh4SyQS6mKs2EiWzg2RJyRIGEEYxRSWj5iEomXCBDacEG6r5q29w9n+eoLM9EVFPZvEq5IOtRNREhRPSZdayk3DM3JAM7RDxDVE1EJcjHp4hwnfcXbtD6M3MT7TmwGwDgl6eN0vowbQY292JtTMRmoKAWdDr3KMdPYc5PGnmbI8ZnOjYTihiH/Yc2FLXYv3/NcdbrVDaDpy89HN07VoQRyI50FOpe+7bHBFSeLa4/W8f47//fIgDAC5+tCq/BwvUrpiaUDNL9DJocDNSAHN7ufQd3R8/OYW6dPQZ0DSZVfVZIBhkm9rnN5kk6B+OJ6QErIWauFCEaKmthxKMpzzGsT2cct2f/4Ng+g7pb3dAU1GQ3ubfBvTpZJYMMCeiqtHDSM1fW4NVZa2CC+txvqw/14S9dfpQmGagNx6l6R3HzGcsLAIQvtIogppKBmae/PicW7NDmA9QGNpEyra6JPFHt2VwHs/J4uOjLa5IGTJtraTBOQ01EDZTV2xuwuTZLGBDHwstpSo6wNcp0fCw9uRat36oRDDMCeWjvTtitfxdniox1W+qDcfpAF0gFNTYa+c2Y3ZsuINRE6gskg8DeoI/NREOea9dH35ucQdQUgmdKFn8bgaXI50PJTM1JumjT58ak7TCX417JyA2bzULdTzVHoU3m+Ws24625jQ+4taHdEgMgGoHcr2sVRg3spk2YKqLj3ritQeoGdV9lIMpV2fTRgAgeA4BT9h2oj0F+UzFVHQtS2Tquw3SrDMZA/obEQHgTJdUX21CXzWvSzu4DukZEf4CqdziJUrVHUgfSRCAZsNADyXjJG3JCQqMLwg0TZwEA1myulWX2sWfJy77NosrLcS4Ns+GizTlH9fYGMR6LZPDr52YCkHOB3HM6j96auxbrttTFLnqcG+6H8ie1TdFEa1TNqCQx81pNVaS61r+9tRAAInYrE/m8feGkmXvLPRk2F0qJW9ODB1w7167DRSSFmsiuqsoaC6gajy3SWiPWFldQ+v6FczKsp75/fPxugSS0tT6LhhwveCMam73dNDQz6JqF4//0Fs67/4PEfRSC9k0MjBcz8JemL5j0ftmwVXBRD76zWDNumt5EYdsMI/sJLyI6P0b266J5GpkTQunWgXBRULpm9wKX1+pTqHPqSGoHynVXlul5fxT2HtQt2pb8vu21eUHZV+WG5lT0pwFugJjMj7z3OTZsrSeumXp/Su9uEkIgqsoKDMhkEeveUUgmoQ3EjrpsPrhPSzZsjRzPSTdJaht4d8F6AMA/3l9itRkoVBh7JNjqmQyICWrABPQ5IM7jmq5fsxlEYl6UnUPvIyrFRqUdimw+71YTkcynrvTkD76zWGufGk6z+Xyws51vDELN6LcZUGKRzYdkSSV7FHspGATDNOASm83Ez1YG51FmAVBpNcT1PjdtJTZsrfcaew/ftQ8AYDz1hLJcccQN10FgmwLtmhjA4EaUkYq+PIrzX7dFN/wEC15Gn8VHyId+1Mg+eOWKo3HcHv0ibohV5dHbrqq8NGN1sE0m9XrxTQrTrXK4lD7oOUHQWUUZMhmGJRu2ASB+87IttSjTFN/mGB9593MAwHmH7oLfnb2v6Iscny1VYYrLpxw4Ve1QBN5O5SqjY3jM1GkHNgN1/ZzjtP1EbvgRfbsEZTa8NXdtoKo7bESfoDxPXvayjN52DXERDRfzaPumHQPEztOjUwW6dSgPLsy3uY0tMImqQFbXCOnHzBRaZdg5Xpu9BtOX12gL4IBuHYLnpGwKQdtONaTFOUH+3SLVhdl8XqaZcC9c5cQeZEoGPiKr2q+w2Azg8SRTbU14Snj4banLGmq8qGTAeXTr0hwhLFTllTGu19DoaNitv2AADybu2rb65rWYdcbs0jNWrVcs2jUxsEoGhl5RiYj0BVU2AyCa5Ou9hYKL3KU3WZDJw5w0Zy0+XVZNxuBeHFTKhw8WbfBONNMT5z+XHo5XrzhKO8eUDBSC/ZllvZdnrgYguB0TaoybZZ4hda3iOsK+VOCc4uifm7YiqOeSDLZTNRZ0LyvVjiIgL0xfhXlrtmgc7fb6HBgjG40nYKZ6EdsQXZzKMuYmNmFjqpwSs28dNhRdq8o1AzGgE7Sde3bCmKG9wuNOyUC3FaifNOfS6hrBmFSYkoGFyQCAHh1Dew5dsPcb3B39u1XFqq6oaykF58B9kxcBAB5653On9LrXToLodJPuzCYx0AzrTjWRPgYqGcxeKeZw3CKpgrgUKsttWUujxnKa6nrFJpHc8NOlm4TNQItY9jBsRKKg12DWNtNWmPe0V+dK7NqvC5oC7ZYYsFAyDZ5IVoqipocG51zj5jPahNINlqYbWxy3FAzBUkW99KIh/8sKhAtG944V2LVfV23i/+GlOQDCdBQKw1XgFnGBTDrG8XuHLou2pFtqcVI6YyB8iaOSQWjglg2G7RheMgrUH/z2SfPBuZ+4mqDPVF+cdK8hOlYbMctzjjLC9VLoHlAsduH9YNEGTQpVni+L5IZF9BZ0qSyPqDTNOoAlbTLCa9WT79nHlM3lowZkeVZ9Tjy3bfVZ6XcfbWTfwT20nE1KTbStPouPl2zC5HnrYu/L2/PX4WPptg3ojNyznwpmg45R2eY45/jqWBHzs/uArjActq2cuSkFlRGvI2W8fXnm6iDOIGzNp/5T49alG/N+PSyl7vBadILlyjRcCrRbYgBEF476XB4yE28AZTNQz+zq8XtoD/R9yR27Fp+4x2YzItnbcVMD345VEZ010yWf/Qb3sNazwaxCd2wSx0WNo3YTGWVtNgzXYl0bkQxCmEbu4X0749R9B5JoYNK+GovleswXr7wsgy8dOFgbz5tz10quT11TuMgO69OZGJDDdmjwIV2c6TWsqN6OXJ64RDru9/w1W6zlSu9N79uVJ+2mPUt172jbnYl7K6A8YEKVmMpiK86zD6qmNht5luqv8t7q06XKaTPI5gwVj+R2N25rCMti7osJW336XnYn0pBSC9Jxh7/1Dqkb7b6DuwOARuRoMLNSE3XtUI5jZXS6a/zq/DjJQEmrqp64V2Gt12aHe0uUGu2bGMhvzoE/vypykt8xaYFuZCoLowwBoB8Rq4HwYbomgck9AHqUcdKXQEwchwjtyHtichW0XMFcaH3jMI91ogsNEWfLmHCHtRMDe1uvzlqtleuqOv0e5Qx1nsqzpLVvGb9tofrDl/YLEg9q4yT2JLXQ//6cfa1xDA25vPW+0WvYXJvFq7NWW92ZCwJHQMD6de2AtZtDKaKzjP1QLXeqLMO5Y4dopzMiGrw1d60wfIZNR7BYSiT/eP/z6FA4x+idewAAfnf2Pk4pWLkV0zFwwKH2SXZf4u4jNaZ/unRTUK4txojOQ7pJkFIFV29vsLqWZjIMW+tz2Fybxeuz18S4f+vqH9W/iX0Hd8eAbh2CNahpZAA72jUxoHrhP786j5SHdYTxiONJGeizqrpOe6BlweIgMHZoLwzqEWaupC+ImhADZG4cdVyc738JfCLodc/O0K6Htm0PlqG/9YXWPwpTTWNXP2yty6FzVZl1IrsMhfdI3fOm7fWRtsOFRJyUzQlVjqqzeN02Mg43cXVFqtKxA8BhI3oHRkTOQzfcjhVlVsnjiSnLsF56mynXVljqqb5c40sCDkHU1Bx77MNwvyhFnDnnyOc5ttXnAhVK0D95KOu21GN7Q84b8KUcDcwNhdQ9+/HjUwEAPTtXOudoQz4aMGZmG/URJEBwzN88ZJdofW6XiFXAn9meJrlZOPMtddkgoHJwT3GP+3erCq6rb1eh7jpj9E4Rqd6cRxTqUl1GboVPl1ZjVQ3NXFz8XCkU7ZIYmNxnnnOctFf/4DidMFXlGeQ58PQnYouFBWu3aByNOQ/Lyxh26tGBHA+5c2UkpRx1hryIvqypPre1eVK1MHf1Zv0cT1vmb9eG8BQRXbTZjjz+weINeG/hBqsKTBWZHOQ5kts9cdSAyBhtm7aoiFcA+kLjIa4uewgde9eqcuzWv6vm9RRGR2c0d1YbnpgijOfH79kvaIM+lzjdeBLQfaG/d8yIoJxy32/MFUGDVHJQ4DC2V1QE1NKXSz+triP0qss44wyyMi6Ensuhp9i27SBntmHbRc68FoWOFXaXaXXeATLDrHn8/UUbUCPTqH/tYEF8DiJpRbrKTL4XHj7MwXzZx79s4zatf/E72v+cyDscvadUBVZKtEtiAOgLPgcwqEcnACLNAn1g/bp1AHiYcoAmSAOAi45Uid7EceFLr8uiahHaJr1wOlWWk8Nq0eFB5K6JCeP3sE5cEx0ryrzHzT7l8OT4xbfSO//9IpGS46nvH4YDd+mJnp0qnP7pqk1z0ppZLC/7wq5On/buHSvQubIsWORshtE8F9lXV9XU4vEpS4P2yzLu66FwEgMy9jrp+kuN05f+82M5DuYkZgrKIE/TWZx4y1sAgEuOHqGpn4oB5yJholpcexOPKPWbIzREKk8e7Vo51zdo8ahofGoKXV8v/e4tU7jBTCUhbQvzVgsm5jtHDotVE+kbSYUjczEvnZySgUCvzlURBkv1rTzsjt6tLxbfeAoG9egYXGt9ljAGEcnA7U2kNA/mOxE3CyiBUZmHlf2o1Gi3xADQ9dcDugvx75Kjh2vEoEwGptA8LnRSDu0tiAgHsH5LHWav2oz3Sf6hDJlxyt+eSgaUw1GeI1eduJs2zq8dPMQrGSh88cDB2n81kQ664VUAwP5Desgx0Uph/0CoElH344AhPfHv7x2GqvIyy2If/b1Kbs7TVYratK+rTtpdk8YottXntNa1Z0D09NvqcpE6egIzN5frVBPJ+6Q24aG5clYTkb2CeJq5nsVFRwyLtKtw95sLyMJrP79v16ogkM8GMRd5EG1M3UmDe8ZDIqA8acK8QGLsr0gXYu08C1T7vzxVz2+kFr59BnUHINRWiW0GAMA5nvxIqLjenLs2Vk2UzYceW3TMtj2ZAaBDZWhMZww4cmQfeZ448Y05ayJqLRXrshNR45oII/ejkd1JVDobt4Y7BFKJVGHXfl20xIL0va+tF30fu0c/fydFon0TAzIF//7eEgDAlw/aWcsDn5EPuK4hzKBJXc/CgB3gv4Z+VvQRviDKC0BzGVUj4MBTHwsVw11v6KmHla7aNdEOHtbLGsimXlilKlAZVdUL0aWKuCbKxn/46CfyOk0vlGi/ES6HA0ulOKxiEcw6LnXAox8s0YLTJs9bF/xW0ko+j0CEB1zeSmosFjWRfKxfGaMvtkqVpyKN7528KBjnWuLmSblBl2Sg5sDzn62K3LNBPTrGGgTNhdME58DrMpgMcBFCjsE9BZNy3J5i4Xjn6mO12BNK5HxjUmqgYLe5oC/xPbhnR+zWv0uQwsNGbxtyeW3vYrXAHTS0FwDg7AMGx6qJTMlA/druIAbhm83Rs1NlJIgym+eRuansPSpamDamhvXz/4gANqtkkIBho5KaLfnk9vqctoGQtkWo/D56t8bv/25D+yYGZB1UhrLKsowWjKSybNKgKC0HCslNtH5rdF/gVdW1WLxetD1Der2s16KZQ9b8NzK/zla5KP71mwfi0OG9g5w+LlSWZ4IIR9v1mVALwZa6rCGZhETQFsDk43rUi2B26QhcLUhNou5xnnPc/LLw+tpjQNegbZp63KcmUs9wtwHRhS3POZ6SdqG6bD5opzNR6dGspGr8n5EAQgB4b6GelZYO44LDh8Yueg1ZOzF44FsHRdozQa+9LtjmVCws/bp2kLEnYgEy00q7xhSmerDPh9qGXNCHeB7RRpas3xYwB2qcnAPvSOIr3DJlmw6bQzbPsb0+amOgDARFoDLjur0hsrMauWjVlhnURd261YZFGcaCpIoA8I9vHwzA7Vpqjlv8jr4H2xty6FgZJZxAmAbbJQ01Fu2SGKgJp57Lo8Qjo9zYG1csFOG5Z+6/k8ah0ECfft2i4uUHi8Xi0JDLB37LZx8QqnMoB6PpcQGctNcAPCp3C/MtopPnrcNny6sj5S9MX4W5q0O/9XF7CeOs1YuDh1IBECUGG7fVY8PWeue+uoqDUZ4Qvzp9L2s9ajBPCioBKO+Y+lyYY0gZeB+6cKzXZ/7ZT8VirzzDKDgH/kW2HDUzVwI6N6ia/8PLc7R2dPuFbuc5dERv7Xmb2F6fw9b6HB5+d3HkWLeO5ZEyAHj8wyVhf4TwXvdfkTwvCOILxiTujfKKueUro52xH4DbbVktUpPmrA1iI1xxBiuqazFL8/UXEuvrs4WRW212D9jnhdoKVaWBVtcBAP98f0n0BIRqIkDESaj8SFtl+gy14xztbo2UoIf20aUIca36wMzNn3p3qbTGLfhg49W21+c02x+NXTjn7ncBAG/NS7OWlhSMhTpCmhK2U0UZulSV49dn7IVnLzs88sAUd6VACUPPTm4rP+fAljqVHyiq5/3nB0uCh/6bM/c2T/e6/yXFJcTzJGw3XDxfmB7mWjeJQW1DHq/NXoParIMYyMXhJ09OAxB6T5j7Idu2jbThzNE7Bb+tqhnyUxFRofaKHA6gpKfL5FajdEy0/n479wjuNyXQFSSJoVr0Nm3TpcGXLj8KgNya0Zg8u/fv6n2Oyq3WvGdylPI8cUzln7JtbUnPNjO+Kib3kr8Lo3hlecYrTb0u4z+iqR7C/0rictkMTCjVq4I4xz0vtpN06UHvsvvlMj3EIcOFyunCw4dhRN/OwXOaI43BSg+vOPuuylGEdHfBAyIbqGmgpSpaZXcb2K2DpppV+5zYLp9GlJug1Rev2yokA01NFGVsyixSWinQbokBACuXq+wB3zx0KPYd3COiVwTClMWAnv5XuXzZNsfJc46rnvwUgP6CqtZ//2LIYX6D+FMH43IYLgtRt2y0qLFciyeNGKWo8XgycB4uDK4JaybGU+jXtUoLxqN+7bZ8QNXbGyLcOw1Es90WVUTjQIDQoKlw+n47BW03kHGUZ/S02YAReAehYph/w3jc+fUDAYQbJgHQNsexPTXVp82ArM771gMfAhAqPgBB0Je4jqgKysyqazKvdTTOwDKmhyLpEULUGioam81AEdMhvTppo6DqGpo23jYGG3E0jbfH7ylcw3952ii8duUxQXtn3PE/ACFzovDB4g0Rnb3Lq5qqakb07YKduneISAYcoUrZBDXWR9ol1X8gJfNN5B2z2SEqLM+iFGjXxKDeweVS+DwtAPICQgRbAUCnqqhIf4/MHR9tP9mDVbXMyabsEUmw1cNhmXM4wlFKnH7729byyfPWYSqJ9DzH8GwyOzSvY3uDbjij0prNA6lbx4qAsKgFp4zmxbd0PdGSfA+QKSiWVeOLUn33dbJ/9b1vh6oJ26Y3W+ui+lu1ANuerM+1VM3HA4b0tJynQ6k0viwJKJVmfAGMjDG8SxIMHrVbXysRMWESg0mz12h2AEDf21hBMQfnHRoyOIyFXjmAYMCC1i1jsG1Ha14jddcGovdrUE+dAQhSTlsuOpquW6j7cnmOf320DCukxxwlBt06VDgNyOsdkoGZ2G7pxui7TKWSYHwpMSg96h36bwpbArade+kRxoBYkC75+0cA7Fz1za/MDX7360aSdiV8rq5Fe/H6aE5+F0yXtN+cuXdwfaZhzUWk1O5YcTC5b4VgHls4SPpybaVpr+VYKIclIkBF+TRpxCXaButLrmwCdDGk+Lf05uog03xHxp6hMSbSgCxtNWOlZwxFUkKvcPvr8wEALxJ1XVxbI/t3QceKMlx5wm6ause0FQTtINybY0ivTuhQUZaIiFDvIwBalKyCLfVKtcw/VGWoPpR3njLW+uwWNu87c76a75x5u1RqdIVhfTrD0BI5z1WL/BOGrYmqifp2rXJ6/NkcSwBga11WU/OpmAnaxtNTV2DJhm2a0Xi0VFWVGu2SGNBIVoquFo7e9g5++4jhkbJNJOmWzzUQCD08AL9bnz4O+8tygVQb7EfUBS6YQWlfGzskGMBfXp9nOcMN5V/uQrnxcj55yaEAQuJq28/A5ioKhASEutz+4NiREU8llVM/Djv37BRbx2z7oKGCW6eSAZ0/t351dKSNVdXbI2U+dYja21h52QDA4xcfgju+doBznnTrUIFZvx4nkgPKSiurt1ttCbR/IPSgU/BJBiry1geGqMR35RNCNXo3eXaMhRLDt2SyQ5/d4j4inSmY765p4zKZODOp4oVHDNPUNHSxNZkmtcib0r0pPdsiigGxV7UN/5DGb7WFrLqknKWNPX7xYvD7kqOitr9SoF0SA4Wxw3Ru7uUrjorUoYvCp9eeCEB4+QAQBmZ5/BrpfwzEEwMNBXKPrhf24iOjBIqCZls9cBexsGWI2oNu4JIEV5ywm/e4aXBUPuU2Y/D/PS3uHfUWUThhVH8rkVAb/lBwuHMfURwyvLf7oIS5mByzez85/rB9ldwQsKcI+HDxRmf7vvE9fNHY4PfBw3vjlH0HRqbJhPF7OM+/Y9IC5zHqXaYQ54IMAHvt5Cb+o+QmOTbJQBl4N9BgK7Dgv1INJnluXx4Tqh7NVPE2bycKk1nq2qFcU9NQD7MenSq1uqqxhet0KbyDJeLfNnyaNdUG5dKq7EBxdkCb1FoKtFtiwBjQtYP+AtsWcaqTVy/8gO4dsPjGU4SB2cKzuYyvNtCMmwDEbliO8fpgcuIm6IL6yEVj8b8Jx4p2LXXPPmCQvzPEq9hcag2bv7cK+KM2nHvPGwMA+MM5+8aORSGX507bCkWfLpXOYwrm+9axokwrz3OOl2aE6hxTZ23iO0cOA0Dvi3t8Y3aJ2gxMfMvgdAG7SjMJ1Hk0oA8ALn54SvDbt3FMj07hpjXmfVfE4NARIQGmNoAKYx8I33OjRNyUDA5NQOApOleWa2qdLRbbj8JctQGUAeUVqCKWTe8kBbWF7EMXjtXKlaNAfTaPj5eEjENjPAYbg3ZLDGyosHjAUN/zpgBVLwFRH2cFm5HvL2QvYturesnRoTj5TWLA61RZHuj0bYv2TV+MLsAq507YRrI8SCYU0TI5OxPHj+qPxTeegh6dKhO5KwLATjIlAhB9J7XEbMY1/2xcyGUHhNC4LV8/RBmVQwP4yuqo3twFxXi41CENlqRtFOZCb404L5JhVOf99U1dDUJ35PMZLW/+8n4A7JKBwk/H7R78fnXWGmdbvid9MknTsIcRONibbJ4DAH9/73NrG8fs3hcVZSzYd1n1N91gyihenKHbcB64QAYAypOVQdkWUQwAe0up6kgjsllpJuqyeZx95ztBuc2QrJA0/1gxSIkBQTHuu7YXsCcRM6+MUadEMx/aX7pJMkBH5ZXP57lmlN5kcfk8bb/w5bGJtICdiNgkJHMi0/2DTfzaEiehoFxOc/k8VlXXYuiEicGxHxj+/wqunEImRLoHu7phm8WTSoEurMputHFrg1FHlwxWbKoNfNaTYHe5eLmIVVxUadSoaSMYOlQ+njio5632kbb37yYG/bqGnLG537cCjeSmC3mw16/DQYLOD6qjLy/LBIFzNri87B68YCzm3XCy6FNSg1kra5yeZjZ8QaoMTUbO5Z2k0oyb6h01p1R8g8JfvznG2bfv3Wos2iUxKKUUZr4j3zhkSLDJCAAMMzhqZYgMG9D/Ul97igUyFP0RyfGYnDXdWlLB5R6qdZ+QmzT91X2cYqVHZaXUDdk8xzl3v6MdcxmlTcmAuimacHnG+AgKJZTq2e3W377PrFoUb3pxtrM9G5SdyUas1tTU4nsyCMyFJNKRuWCftb9f3TecSKFdqsrRt4t7cTXxe6K+U1OBpmsHdO6cLoQj+ob3tiJww41Xn5mLqbolFx/lt5dR7z8Kdbsu+6f/3rtw/ChBFI4OdvaLjv6jz+3GYyC0yZiEq4vFkUXBF8DWWMS7CLQjFGT4TdgG/X/X1w/AeCLqArpk8Jev7u/kzpR72ufrt+Hqp6ZFAtP2t7ibmQu4DUn1zFRnHCeq2vLYmO1c/9+ZWLZR97ZxeROZ5onrz4jnjsy18/g/veWs27kqej0uTthGA6/2GHNNqPwy67bUYXcIDvma/3yGt+eLxHxUnUJBib9rcTNxkMXdlYLq8SvKmMbJ0+03v2eJXKcLerArlyEZZEmwGE01PfGzkAtXNgO6mVBSqIVxcE///ahzeFblOUd9Lq8FBsaBqkeryssw+adfCCQUuqWowhfvetfZls1l1sQJo/prLtU+201j0S4lAyC6CP7i1FFWVYrNUEeRNaIjTfUJ5c6t3gdk8ozbe4C3L4VHP1iKU27Tg796mh4QSDZxzDXPxZ3/lbjVxRmYX58T6oTfuOoYvHbl0eGY5MtP3ScVXDQkiZoo2AfBcck+jkrbKS3mltmIpy1i3AXlNriRpLGoJ3OoqtxOaOk9oLYgfWwh/vbNA7FzL78L7TAiGVSUZTTDLiUM1KaisN2SIG51Ta0moT5G8ia5bESKAXC5HCeBLT8RDXp866dfsJ6nVK8UN31xH29ffzeyC+wsYzWA0AX1yic+1dRbLhxhZke1wJxbLoapFGi3xMDE6Q6O/IfHjfSeZyaXMx8wlQyqLIFAlJiYkY8Uwx2GZQXTTRZIplow8d8fHBFbh26ookDz3VP969A+nTUu0pdXxcWNU4J57WmjrHWU15G6h0ntDCa+61hoFWxD7OwQ6wdYEheq59S7c6iSmUwSj22ts9sh6Dw5c7SdGNOxucZ0+fFiPleWZTTf+4qyTLBxCxAf5aoeI+1z7uotqN7egPlrNiNrcNxxHlxJAt9cMN9BQGykpOCyl9ncbL9y0BBLzRBDfARWqolU8KLPTgXYvQ5Nw7ip6j14uF/aawxSYiDh4qLjXgqT4+loeNnQB25Tn1BuzOc/bEspHQfTw8KGYhiNiB829AVN5YmxwSetmFGiCjQIhxIWCsUxlZdlUJZhgZvqX99coHFpNhULNbR36+DfUrCQ+/Xi5UdGypSxuj6XB+cci9ZtNXIJ2Tug88Q1F6jUMsCxQYtS13zvmBEao0LVRPNWb9ZiKGzYqbu4jz8izJJihI7/01vY9ecvaPVNN24Fdemhy254TNnJ4mBL39E7gftwMfBtOUlTXQPAqF++5G3LFlx2+9f21/6bDKQv3qOxaJc2A+sm8Y5FKi5mII4DpXr7zbVRj584F8sk/Tz1/cPsfSdRE5EF5NZzRzvrnXPg4MDNdqTFuErzMf3Z046LuB4wpEeQedIEvfajHBt70GutKs+gLpvDR59vxO9e0A29NkJlI24UNB13IWK6zRCoNg36+PONuPnlOUEqDYVejrFQycD5XEnxLg4OVqmkTIJSUZZBQy6PP7w02xu0pjC0T2dM/ukXtLQjvSwSIwD83yl7av/F8xGER7ldKonxw0UbsGjtVpx36C5YTmxKNIutwln7D8J/Pllu5dY7x8R9FAvfOyXSViSXbGwqWZPZoTYKG3NRSrRLYmBDsZJBvDoiPL7FogKw+YsX2o+NMwKSGcTp2qYiSW04a/9BATGw9UdTefi8IVz3+V+XHOZREzmbC5AxiME9kxfhnsnRiGabzj9ug/HdieheCDGwzR21XeWtr9nTf5xuWfSAkGkY3qdzopxHLucBJTGZ6oeKsoyWwlyBbuNpwrRJuObbMEPF2bNTZZDbSBlf1bxQ3nK3GJLJn748OtLuLV8ZjV+eOgo9LUSo0Cjdrh3K8d7Vx1mPXXL0CLFlKfwutgyu9OP2QEJTWrrk6BGR9qmzxh4D3O9nKRC7WjDG7meMrWGMTSdlvRhjrzDG5snvnrKcMcZuY4zNZ4xNY4wdQM45X9afxxg7n5QfyBj7TJ5zGys0u1exkL2oF9a16NsC0SjUIr3/kB5BugoKmtVytCV/0A1n+Q1WCl8oYt/TJO8DreIjfIpoDe/b2aqPVsZvX5oEwK5m+fgXJ3hfXl9enBNGCU6/P9HPu4ywADDQoj4pyzAsvvEULL7xFOs5lADEpT3Qj4mDNH4izt3X9QxUqdogyd6ft2kAwHypelF7Jyi4JOBvFmAcd/V/kGHPOu8w0Sa108UxRa75YSMExeD9a45z2ll8e1JT+O7/pcSGYcP3jhlhTfFiqp2bEkkkgwcB3A7gYVI2AcBrnPMbGWMT5P+fARgPYKT8HAzgLgAHM8Z6AbgWwBgIVvkjxtiznPONss53ALwP4HkA4wDoCscmxNcPHoKH3/3cydXEcRjq5d2pe0crh6n2HwaAXXpHjcB9ulTh2D36oX83v37/6vF7WBN23Xe+O0AlCV2ldXxCjroPtmR+gMjdM+NXJzlfKAVqpASEF5dLvaDg2wD8zq8fgI8+36gRWltGTYULPZyuiQOG9MDHSzZpRLVQ241JYOLOdzEfX9ijH3543EhcdLh7/Emkli8eMAhvzV0bBE4F/brmfwG8mSta37TDfO/oEdhnUHccTjzvirGJxeG2r+4f63gR9O8h0i5u34QtnblCHDNn89gCmk7dZUNsT5zztxhjQ43iMwAcI38/BOANCGJwBoCHuVDKv8cY68EYGyjrvsI53wAAjLFXAIxjjL0BoBvn/D1Z/jCAM7EDicF1p+2Fq8fvWXSO8EOH98ZVJ+7mdC88Zd+BuFJuauPC/XJ/Wx/K5d7MG7bWY+ywXti1Xxecd+guiUTHCz0LCL1u3z1QOmtfCus4QgDEZ5u0wUfUKsoyiRLPAcDEHx5R0HNWQ6X9mxz0dwvMIOm63jG79MSZ+w9yMh9lGRabHNCnnlM4Y/QgHDq8d2SLVjODqULSmIZCwBjDkSN1249PmnvzJ8cU1Y/LQ9AGX0yOK1+YCbW1polfnGr3gEuCYtO+FINiyU5/zrnyH1wFQFnlBgGgSb+XyTJf+TJLuRWMsYsBXAwAQ4b4XcB8oMtRJsMaJYplMgyXHet2P+1QIQJTbG6lheLs/Qfh3rcX4aYv7hvRw7rgUn1QPHLRWFz/35nYqYfdAwUIvVwKScJng7nBjs+dluKDa45LFETnQ6GeGMqDhxIQqsO94oTdYl2PTVSV2efakN6dCopXcGHeDeNj7Vi2vbrXbNbjMO7+xgEYt/fASL1CcdtX94+vBD9TkITJKAZ9ulQliujt160Dnr70cIzsZ/dki4PP7hKHHaU1B0pgQOacc8bYDsmzxzn/G4C/AcCYMWOaKbdf4YgL/kmKCeP3wBcPHJyYECTFkSP74pUrjvbWOWBIT7xx1TF+P+sEMO0mSf2mbQtYUnzpwMFatG1SXHni7rj8sanaIkBfTpsNKA7KgExx4C49NRfNxqCiLINS5DJT6TMai6TcuU9is8W1lAJPfPcQHHvzm4nqFvOsAZHmvrGYef1Jib0OG4NiicFqxthAzvlKqQZSoXzLAVBry2BZthyhWkmVvyHLB1vqNzkKobcdKjL48phkRqSmRHlZBnt6PH6aGq6MqoWAGnrn/mZ8k+iKTfzhS/sVdd4Xdu9ndQpQKOZZ2HTT//6e3TW4OVEMR3ruQTsHSdkA925rO2IsSTC8bxfc8pX9EqsZ4zDvhvEYacRX7Du4h/ecdyYcG5t+Pi49eqlQbC/PAjgfwI3y+xlSfhlj7DEIA3K1JBgvAfit8joCcCKAqznnGxhjNYyxQyAMyOcB+EuRY2oyzP71+OYeQpvEjiAETQkRYZo8uRug66Z7dqrwBui1Nvzu7H3QtUM5fnzCbtiwtb7oKHCKpPr6YnHW/o69uosANcI/94MjnIF/FDs5todtDsTeacbYoxBcfR/G2DIIr6AbATzBGLsIwOcAviyrPw/gZADzAWwDcAEAyEX/1wA+lPWuV8ZkAN+H8FjqCGE43mHG4xRtE5d+YQTumLQAPzlpd2cW2MbguR8cgVtemdtoldnvzt6nJHr5UmBwz45B4kC1RWmhYIzh56cIY2mpuNn3rrH7/rdUTPyhSOfSlJHCTYUk3kRfdRyKPCXpRXSpo537AdxvKZ8CoOmSdKdocfjq2J2tic5KhZ+ctAeuOnH3JlMv7D2oO+5L4AEWh1KpJ0qBX5+5d7Cfdly20x2FJM4PLQ2tkQgotM8I5FZjem6b+N3ZybeyLBY70gujWMSlwdihaAHvxE/H7Y7fvzgHn113YmRzoRRNj/ZJDNA6FosUbRPD+nTGonXJc+jvCKicOjsl0HM3Fb5/zK74/jEiUteV2C5F06HdEoMUKZoLk646prmHEMH+O/dEl6ryxHEBKdoeUmKQIkUK9Oxciem/Oqm5h5GiGdG6fftSpEiRIkVJ0C6JQQuwlaVIkSJFi0K7JAZAYRHIKVKkSNHW0W6JQYoUKVKkCJESgxQpUqRIkRKDNon37gKu6w5s39TcI0mxoyDTbadIUSzaJTHgSTbVbU3IZYEta8L/Hz0ovjevtFZv9cjn29fit/ht8YxdWDQZuL4nsPRDd52alcDyj0o/thRtBu2SGLQ5vPxz4I8jgdoa8Z/Jx8pb+YI56znguh5AvRGt+5f9gZsavxFMq8Dn7wAPngK89Qd3nfmviu/Fk911/rwPcM+xpRnTtg3xdVK0OrRbYtCmslHMfFZ8120W340lBis/jS7ALtSsAN78AxAnbU15IJRYbFg2Bfj9cH2hee16AByoNvbW3bgYqKtJNj4Ts/4ruORSItcgPk2BmhXie92csGz2ROA3A5I/IwDIl2h8n78D/H6YGEOKNoV2SwxisXEx8P5fm3sUyaAW/Yzc4kpRunwRmUHrtwF/PQp48lvJ6j95ATDpN8DqGf56z10O/PdH7uOvXQ9sWy8WG4WsSKmMsgR5atbMBj77l79OrgF4/BvAw6fHt0fBOTD1n+Le2PDnfYHfkVTZT18K3GXscPXuHX41ThwYeVVf/ZW4N5uWqAEW326hWPGJ+F7kkUKaE7kGYMr9xc39do6UGLjw0OnACz8FaqubdxzbNwLVMZu/bVkV/uYcyMjFM0/0zLmGZAblrWvF97yXxfe7dwhjdK2DE2+QC2RjOc9FcvtBKrKpxS7u+gHgzoOBf1/kr5OV+91uXOyu85cxwIf36mWLJwNPfw946Rr7OZtXhIQLAKb+HVg9Xa/z0jXAfcf7x2dDIN1RUZbrZUoqM8XdvxwIPP8Td9t1W9zP1QUmGQ4uF9t8Hpj+lN+G88hZwH++Zz/WUFta+897dwLP/dgvha6dA6ydm6y9VZ+VZFjINYhrbcFol8QgER9Vu0l8NxWHUb8NuPMwYOkHevnW9eGiBQC3jgZuGZWszZf/D/hVj1BCoKqLX/fR9ezLPhKL/PoFehvryEuyYFK4AK6fF5bn8yFnqPqiL/SUBwQRKQY21db2Eumo1X3N1YdllDBwLq5z4pX6eXVbxPeOMMjXb7UvjnShV88oKDOIg8L6+cAHf3P39ceRwI0xm/8smyLmyTr5/IPnLRmNB8YD/7oAeO06dxsLXgc+/We0PJ8HbugvmK7GYNVn4fi2b9S/bbhjLHBHgv0opj0J3H1EqIYFgGw9MPOZeLWoib8eJa61BaP9EYMta/Dm1jNx2hpDBbRospj0K6aK/0osp8Tgs3+FCwMArJxW/DhWTwfWzIhym38YDtx5SPhfEaUk+OxJ8b1MqiPoomfi00fF9/zX9HK6GD9yZvi7nKQ2nvhj4KFTgTkvABmZ65CT+/Tc5W4uOg5PXgA8dJpetvh/xbVlYtWn+v8Z/wFu3S80wAZqFwMr5Xlr59iPlwoNtcBvd9IXVrXo2BafO8YKG4iC73mb7QGhVOfDvXIPq9sPEu/H81fJduQ8Wfqe+H7vrvi2TCiCMuU+fz3OgamPAg3b7cfvPgK4fUxYl34DwCf/AJ44v/DxrZwqvjcuCsve+j3wxHnA3JfCspd+LtSIPqyZWXj/NqyaDmxeXZq2DLQ/YvDHkQCA4zc8Crx7Z1g+S1L/6f8WE0lxFmohWP6xUENMvEL8/+xfwF+PFAsKxZpZwMeP2Ptev0C8ULOfD8vopFWG0g0Li7gwC2xGTdXfh/fIOnX68R6Gl043uUdshujtlQj+9PdCYlAqA2q+AVj0lq6+WNUIokvRiewsdu/xwIdyEVKqgBqHOkoRTroo2HBdd7suPakaRDkA/O9WQRgaakNO3OUM8Pg3gNVyoZl0Q3wf2br4OlYYxMiUmJMQIhOKGMQ5Oix4DXj6EuA/3wUmXiUW33mv2Ov+78/ieytxtX7m+8DMp5ON6bruoWS4drb4Xj8/PD7lAfFNF/d3bxfvQhJMvEr04ULNSrH2ZB338+7DgZt3S9ZXgWh/xIDipavD30qcfuc2oWpRqJELtNIBT3tcLNbKYEonCiC4+mcvs/e3cJL4fozuJEqJAVmMbtgpyRX4YdPjK25HYfMq/f8jZ+n/lYRke9m3bwSWvOs+3hhMezz8bd5jBZs3jal20+oTTnjZh6ErplrYpv/bcR6RBs37ZeKhU6NlSe0pswmXf0N/8VG2FB9nuWFBtIwSICr9mcQfABa+mWx8FC4u3Ydf99O9kPKe2AkKFUMz8xnBxLx7O/CPc/znrJvnP26Dmh8f3gvctn8oMU4nDN+2deLbfI+SQjFh056wS3t/2gO4aSjwm76CaOxAQ3j7JgZJ8NFD4vvZH4RlT38fyEpj0LICAnmWTQl/33eC+NYeNhXhtwpdq8KS95L3o2Dj1v92TMiBRvqHMIZSKL20bRGhsL3YisOi+lvatw+f/D38rYzagD7e31oIprqvNrjUImqx2bY+LKOEgwZ83by7u/1Iu3LcVFrwcea0f4WPHxbfiku1wWYQbyCE8u9nkz6k/YVeU6HeVUDyxZAueLk64LGvhcFvnydU/yV1of3gnvB3ZefocdNF2TToUoaGSuf1ljmr2ip2sX7qO8D94+LrvXlT+LuJMwq0P2Kw88GF1V/xcbRsybvCawEA5r5gP48u5ApK3UCRJRPSFJeXvB/+LsYg61r8qATyfoyud9Pn4tsltirYFrkP7xWLDrWz/G6wvx2FgfvZy80XwvaCuMbquh9K/TOHPMvfDgx/+zh732LwzPfF99t/Cst+089dPymhpOjjIE4uA6qScsxrKlR9ZJPWbOpNW7tK+qSeej7VCZXMfFD2DMDOnNyyl3BBVlhkSESue7DvV6Jlyz8S7sNbiP5+09Jk41RYajB4trk098Xwt88oXgK0P2Jw0cuYy4Ynrz/ksOL6eeSsZLpiyvGZnOEupO9Zz6JgON0GTf1vgnHGqYEUUTN98dfPK06lsMvh9nJTQqnfGl38XYuqaxw9h4nvrMP1z2do9S2iilB1S6jyK4YY9N/LXl6zwl6+Vao5qJQK6F5kxWLyn6Jl6y3qGkUETOcFF3z3xfVMlQRkPh9qV/nnl/VjVE1H4XrGq6fr/f95b/c4k8B2LSuJ08PCNxrXfgzaHzEAsBsnHIztAVBD4xIZBFXl4VxcqHZ4p7hgSgauxSkpXpA+5lvW6uXmwl5HOLSDL7G3Fcc5qrGqoCSFTDmwfEq0vgtdBsgxOvozVV/5hijH5CJcrgUyzjDsg8voDADLpP1CeXnFoVCff8Ctdzefg4LyqjHjIFiCpWD3U/zHbQR8poWJOeF68b0gITFwPbdsvdvDS9178/32MVWv/8Ze7jM+J4lDSupAMO0x//HnLk/WTpFol8RAw8bPo2W2F4wumBWd7G2ZBqFb9wP+/e3kY3njJv2/yb0VC9MOYHLS9IV5/257G3E2g4AYGGq1bB3wzKXh/zg1nQqgcxGfjx4wxpWNjs1FRF+9Nvy9D+EKVYBdEvQaof+nLoZJ0HcP97HuRIV23LX2OubCYs5VZQtwzZ1dDhXfNLEhAJRV2ut3I2PqM5LUr4rWtUWK23T3CjYbiQ1UvUoJ0poZ8YTcJ9EOP0b/X6gaZuiRSBS11JDQ5rGhEUxJCZASgw7domV0cVSukxS7j9f/qxd08h+jdZNyhYCY3BRv/T78Pexo+zlJuMmthmSQNbilJLl6YiUDeVzFOCjkG4CeQ8P/GcuCYQPVE1d2DX+/fUu0/cj1JdB/J/HwsXl7mIuLy7YByMXCbNPDJdJjgw6MHs/ndZdJwEIMYqKs1QKoPGWC81zSFDG60piXPU+LVLUSYVsEr60vFdkcB5qjae1c4J2/+Osvt9j8FDr2Stbn4LGO83voRNckLgr0Hf26J2VKxx7JxtNESImB7eWkk3XIodHjpneDWlhcYqaCjZvqu6f/nGBMjsXLJzqO/rr47tBDLzeDq2zeEibUAusymCoCOvMZvTyX1Q2LNqmLXpsiem/dTI57uLtcQzRiOIl6zYwPscHWjklobE4BCl0sEae+a3nntvA3DfKj55rE1jTaqvG5VHPq/puxG0niRGiKh3w2ykTYrm06Wfyu3eTui+eSeeZQw/VOo91pudV7Ve0x6vqexYB9w99KujEl6vptevS0+Z4p0Dq+9Ba2bLA243UTISUG+Zz/Rcg1RHWW1MKv6iSBTdXCE7qmuSYu9Y2nHDQQvlwqilSBqm2AcBKbqgOtfzn26x3clGsBpveq+87267idpAZQahRKoHL17vD/fENUP27rg7pSUpsQAIxx5DRq2A6sI4vPwNHRZzj1H/Zz1dgi40g4V2w6/FydCDKjML164gihy8ZgGxdlenrvGm1n6ft6Wdy1MSYkw1y9fX+GQgMXfSkhlERPPYxM+OIcKLFU8+mTh/U61LGgc1/3+Gc/F/6mi7vJHL5HgmC/Il2r0ziDpsUXOj4R/qlZIfL2uFC9FPhTDPeeb4iqL0xMdXCQ6uV9+lL7cQDYbZx9gTNfhvE3hr87dBcBPj6XPYV8VuQ1+uNId50kBmSbPzhVnW1da3fVpXrfE3+tH6voDICLl+Ll/wvL+8goTBpF/uWHw7GYoJzy5RbuzBY5vHkVcDtR1ww7UucOqZtihaEb7zFELHi/N2wMLqJOF5LdxgEVFskgSWR6tla/T5Hj9Xp21+OvE982wkX3P7j4zfD+AmLOPGmkePBx2qMlESurFHPp172jdWznUyn2wpeBg74T/ldBnDbkG/zpYnbaX7/nPrdQNfdNCVQFXAKS0Yl5R464Aug+CDh6gvjfi3g1mobwPU8Tc1wRrMakvkmIdkkMsqwCMzvL/D8PeAI/+u/t9xbpJInITUOBV68Ly01O87ruIpxe4VvPA9dVAwecLyb7dd1FpkuF776ln19WEeU6PrxXj5Q+/XZgf8I11lZH1T9XGu6DZ8vsnPceF9W9mmN4cYKfsLzzFz0I7PTbo3XUIr11nUjN8e9vA/8mL/eAfYFyg1s6Uqb/+HVvfYzH/kJ8UzVEt0GyH8tLef9J4W/TqDnlPj1y+AwZ03GXoSIsqxIvfD4nApzuJMbwnxvS46YlwJyJYcQqIDy1XAsmZUi+9jjQsWe0zt+OCX8f/TP92Kl/Ft/3Hq/fJ9Nu8cZv9eyuSh9ujmvOC7rbc1UXYNQZYt4CUekYiEpodL4c+3Px3bBV54AB4NDL7GPYul5syqMw5GDgxN+EHm8vTgiPnWAwEWtni3QxLpR3CPtbN093C/3eu3rd2k3RuW/adMqrovPujRv18/Y7V3x/4Wpha8o3iPd09kSRZ0phgiSA6+YKT6bruuvXcpFh7ykRLNbR9oFRWy0RvafeItLfKpjud995XeeW9v+6yCNj4qQb/Mm3hkoXvI8fsh+nRskjfhxKHb7FePTXwjE21AIPnhyt08lQ8XQdYG+rU2+/YRQIFwXXuIYeof8/egLwppRc/jAiWh8QYzfh8rAwVT1AqGdv2CYSp9HFwsTZ94goUBtsqQz+b61ICc3zblWZD6O/Lvzla6sFMfn0UWDS73QDLQU1tPcbFU1H8YVr9OhU9Syph85RPwWOuNweqa2giNXDZwD7nisWpg7d9XEpwuvChCXAjUOEytKV0sMVa9FrhJCiAKGm3LAQ2OkAu+QACInJZu859DLgFc84v/cOcBeJ21FcvW3u9h8lVIIqytpUxanrpSirFMzOH3aNOjQoKGkWCOMHzHYAcf992DlBxtUi0C4lAwB4YoBFlzjmQvHSH/J94DKLYWqnA/T/n1km/pl3ARUdgZ9Y8sUAoRHNhZ9L18ofThWeB0qM92GPU8OEZoMODIkNxS/W6a5/X7zPLdZ+X+qC/28t8Msi0kf/cGp0Qn/haiGa+2C6Jo65SJeYFL73TjQq9ZqVoRH7ifPchOCX0n1w3y8LCc3ETxbYXYfLK91xI+dJ3/W9vyi+x90UrXPmnaF94fpewm5jIwRnyRxZFR3DMnM/h859xbdSvXz7Nbu959ifCynosinAjy2bD427STcCT3tMSJPmuA73bEoExC9ePvzw43BLz7sOFelEbITgHOJSfIChnrp8OpDJACf9VnxMfOlBoDdRgV7lyHUFhO/fd9/UGR4K83p/ukhEMy99300IgGTbK8Y5lNgYphKhXRIDzoH3exic808lB1peCYz7HdBn1/DlBoD/WyMe5vkySvGKWUAnQ5QfuF/IoXfuI/SIFF138k+Ib00MF4Few4CRnjw7FF+xLJh08l9XHS6011WLzz7nRHO1AGLR6CIXm/LKkMhQHPaDaJnCXmeLsVMpRL1UrkAoQCc6339PvNSn/kl3SwWEvr//XsCQQ/Tyyk5AF0+qB4UMmfIm0fzWRPHc+hkv5ASPPvnQy4Dh0gPq7HuAi98ADjEC9woR6/eTBsYO3YQ0dfEbwIhj9DpXScnlzDvEvR08Jmpwpv74fUbqMQxqTIdcAoy92D+evc6OEmmqIrnGERBGoQiwCcUYJSEme5P8Skcaunt1bYdeKj70vQWAUWeKuXz+c0I1qua3ieHH6ETYhi9Kif9nMj7pwpeiErcN5j2wueYCgtFx4aJX7S7HJUK7JAYAwOki95OF9gd6rDTE7TY+1GUPO0q8gN12Ai6QuWzGXAR846monn0cMehe8jbwY0PtdLHMjdJ7JHD1sqhqReFCR2DTcb8ErpxjJzCXSIPoaRY1lgLlsE7+o5AezEUD0PX/11ULvS2F4oQ79gLOuV+vS7mrb1siTg84Xyxu9Hn021O81ID+clxXHaoU6AKiCElPI/22icunR8uuqxb69ss/C+8/fVG/+1YYi2LaXA78llAJKmTKQunnumox9quXhWL9mTF5oExOVElT5zwA7Hq8WIgu+8j+vPcjmXB/PEN/DgqnSdfV854Jx5TxLAFn3wN86YFo+bdfEwzIj2f4g8oAscDSPpQR+ov3hddxVsz2spcZbrIVHcJo52s3Re/HOfcLl/AvPyzuqTo+7Ehg3y+J39SWMnC0uCe2GAA1Z/Y8PWSiABETcF11lCmx4UfTovf5y0aa+599Ltqj9ah24fTbm0w9FIBz3io/Bx54IC8Wh9/4Gv/x459wns+Ljw8znuZ8e3XRffFP/sn55jXFn6/wzGWcX9ut8e2YWDqF88l/Kn27Ljz6Nc7fu7uwczYs5rx2c7K6D50u7tOiyZzPe5XzbRuSn0tRs4rzRW9Hy+e8xPmdh3OeyxXeZi4rxnZtN86n/4fz2pri2nFh3XzxKRT/+b4Y030ncZ5t4PydOziv21pYG7kc5/eP57x6BeeTbix8HA+eFt6ba7tx/u5dnK+eVVgbheLz98T1lgo3j+L88W9yft84Mf/iUL2C84Zaf50lH3C+fmFpxsc5BzCFO9ZUxgvdvq2JwBgbB+BWAGUA7uWc3+irP2bMGD5lSnHpGo646XWMHdYLf/ry6KLOT9GCkc8LQ16lI2VICjvWzIqqx1K0OTDGPuKcj7EdaxFqIsZYGYA7AIwHMArAVxljCTf+TZGCIJNJCUExSAlBu0eLIAYAxgKYzzlfyDmvB/AYgDOaoqPT/vI2VlU3MhtoihQpUrQxtJQ4g0EAqMvGMgCR9JaMsYsBXAwAQ4ZY/HMTYETfzti5V0d88YCEm6ykSJEiRTtASyEGicA5/xuAvwHCZlBMG38+N8bXPUWKFCnaIVqKmmg5AOqUP1iWpUiRIkWKHYCWQgw+BDCSMTaMMVYJ4FwARezzmCJFihQpikGLUBNxzrOMscsAvAThWno/59wSP58iRYoUKZoCLYIYAADn/HkAlmQxKVKkSJGiqdFS1EQpUqRIkaIZkRKDFClSpEiREoMUKVKkSJESgxQpUqRIAbScRHWFgjG2FsDnCar2AbAutlbrRHptrRPptbVOtIVr24Vzbt3QodUSg6RgjE1xZelr7UivrXUivbbWibZ8bUCqJkqRIkWKFEiJQYoUKVKkQPsgBn9r7gE0IdJra51Ir611oi1fW9u3GaRIkSJFini0B8kgRYoUKVLEICUGKVKkSJGidRIDxtj9jLE1jLHppKwXY+wVxtg8+d1TljPG2G2MsfmMsWmMsQPIOefL+vMYY+c3x7WYcFzbdYyx5YyxqfJzMjl2tby2OYyxk0j5OFk2nzE2YUdfhw2MsZ0ZY5MYYzMZYzMYYz+S5a362Xmuq608tw6MsQ8YY5/K6/uVLB/GGHtfjvVxmX4ejLEq+X++PD6UtGW97uaA57oeZIwtIs9ttCxvFfOxaHDOW90HwFEADgAwnZT9HsAE+XsCgJvk75MBvACAATgEwPuyvBeAhfK7p/zds4Ve23UArrLUHQXgUwBVAIYBWACRArxM/h4OoFLWGdUCrm0ggAPk764A5spraNXPznNdbeW5MQBd5O8KAO/L5/EEgHNl+d0Avid/fx/A3fL3uQAe9113C7yuBwGcY6nfKuZjsZ9WKRlwzt8CsMEoPgPAQ/L3QwDOJOUPc4H3APRgjA0EcBKAVzjnGzjnGwG8AmBckw8+Bo5rc+EMAI9xzus454sAzAcwVn7mc84Xcs7rATwm6zYrOOcrOecfy9+bAcyC2P+6VT87z3W50NqeG+ecb5F/K+SHAzgWwL9kufnc1PP8F4DjGGMM7utuFniuy4VWMR+LRaskBg7055yvlL9XAegvfw8CsJTUWybLXOUtFZdJ0fR+pUZBK742qTrYH4IbazPPzrguoI08N8ZYGWNsKoA1EIvdAgCbOOdZWYWONbgOebwaQG+0wOszr4tzrp7bDfK53cIYq5Jlre65FYK2RAwCcCG7tSWf2bsAjAAwGsBKADc362gaCcZYFwD/BnA557yGHmvNz85yXW3muXHOc5zz0RD7k48FsEfzjqg0MK+LMbY3gKshru8gCNXPz5pvhDsObYkYrJYiG+T3Glm+HMDOpN5gWeYqb3HgnK+WkzYP4B6EonWruzbGWAXEgvkPzvlTsrjVPzvbdbWl56bAOd8EYBKAQyHUJGq3RDrW4Drk8e4A1qMFXx+5rnFS7cc553UAHkAbeG5J0JaIwbMAlBX/fADPkPLzpCfAIQCqpUriJQAnMsZ6SvH9RFnW4qAWSomzAChPo2cBnCu9N4YBGAngAwAfAhgpvT0qIYx4z+7IMdsg9cb3AZjFOf8TOdSqn53rutrQc+vLGOshf3cEcAKEXWQSgHNkNfO5qed5DoDXpcTnuu5mgeO6ZhPGhEHYQehza/HzsWg0twW7mA+ARyHE7gYI/dxFEDrJ1wDMA/AqgF489Bi4A0LH+RmAMaSdCyGMWPMBXNDc1+W5tkfk2KdBTMiBpP7P5bXNATCelJ8M4dWyAMDPm/u65JiOgFABTQMwVX5Obu3PznNdbeW57QvgE3kd0wH8UpYPh1jM5wN4EkCVLO8g/8+Xx4fHXXcLu67X5XObDuDvCD2OWsV8LPaTpqNIkSJFihRtSk2UIkWKFCmKREoMUqRIkSJFSgxSpEiRIkVKDFKkSJEiBVJikCJFihYKxtiXZAK5PGPMufcwsyR3lOWuBIg9GWP/kRHGH8hAs2LbKjh5HWPsQMbYZ/Kc26QLa1F9eO7Jj+W9m84Ye5Qx1iHunJQYpEiRotnBGDuGMfagUTwdwNkA3oo5/UHYcwFNAPAa53wkhOuyygJ7DYCpnPN9AZwH4NZGtDUeIl5iJICLIaLOwRjrBeBaAAdDBK1dy8J0JHcB+A45b1wxfbjAGBsE4IcQrq97QyRAPNd3DpASgxQpUrRQcM5ncc7nJKjnSu7oSoA4CiKWAJzz2QCGMsb6F9lWQcnr5LFunPP3uPDrfxj2BH9J+gBj7CeMsQ+lxPArMt5yAB1lBHgnACss16QhJQYpUqRoq3AlQPwUQuIAY2wsgF0gUkgU01ahyesGyd9mecF9MMZOhJAWxkLkvzqQMXYU53w5gD8CWAIRwFrNOX855vpSYpAiRYrmAxOb30wFcC+A01m4oUxJN76RXLiKsL0RgrueCuAHEFHIuSLbahIk7ONE+fkEwMcQyfVGSnXUGRB7RuwEoDNj7BtxfZbHVUiRIkWKpgLn/GBA2AwAfItz/q0SNr+aMTaQc76SkQSIXGSUvUD2ywAsgtiQpuC24E9ed4xR/oYsH2ypX0wfDMDvOOd/pQNljH0JwCLO+Vr5/ykAh0Gk1nAilQxSpEjRVmFNgMgY6yETAQLAtwG8xY1U6knbQoHJ6+SxGsbYIZIQnQd7gr+kfVzIROp0MMYGMcb6QaiHDmGMdZJ9HAeRWNCP5k6OlH7ST/pJPxBc9ING2VkQ+vE6AKvlYgoI1cfzpF4kuaMsdyVAPBQiGeAcAE+BbFFZRFsFJ68DMAbCU2oBgNuBIEdcMX38SJZ9BuBdACNk+a8AzJb9PAKZRND3SRPVpUiRIkWKVE2UIkWKFClSYpAiRYoUKZASgxQpUqRIgZQYpEiRIkUKpMQgRYoUKVIgJQYpUqRIkQIpMUiRIkWKFAD+H0QCPOho5eFZAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "lc = sim.simulate(var.counts, s_ir)\n", "\n", "plt.plot(lc.time, lc.counts)\n", "plt.plot(var.time, var.counts)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Analysis" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compute crossspectrum." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "cross = Crossspectrum(lc, var)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Rebin the crosss-spectrum for ease of visualization." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "cross = cross.rebin(0.0050)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Calculate time lag." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "lag = cross.time_lag()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Plot lag." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZUAAAEWCAYAAACufwpNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAvRklEQVR4nO3dd5xU9b3/8deHIghSlF1EKYJKURRXXRsBgoEoGogNscWI8RfUVBJNYkmi3tjuvXg1ibGgiVgjBKIRYoWIKEUFXRHLIkUFRLqAgkj5/P74nnVnlu07s2dm9v18POYxM6fNZ87OzmfOt5q7IyIikgqN4g5ARERyh5KKiIikjJKKiIikjJKKiIikjJKKiIikjJKKiIikjJKKZAUz29fMZpjZZjO7Le54RKR8SipSJ2b2oZkNroeXGgWsBVq7+xX18HoZxcw6mdkkM1trZhvNbIGZjUzza44zsxvT+RqSe5rEHYBINR0AvOsV9NY1sybuvqOeY6pPDwNvEc7DNuBwoEOcATWAcy61oCsVSQsz29vMppjZGjPbED3ulLC+W0Jx1lQz+4uZPVLBscYBFwG/NrPPzWywmV1vZhPN7BEz2wSMNLM2ZvZXM1tpZivM7EYzaxwdo7GZjYl+6S8xsx+bmZtZk2h90hVXdPxHEp4fb2azzOwzM3vLzAYmrJtuZn8ws5nR+3nezPIS1vdL2HeZmY00s2PMbFVJfNF2Z5rZWxWc0mOAce7+hbvvcPc33f2ZaL+u0XsZZWafRO//yoTjNjKzq8xssZmtM7MJZrZPFfGNAi5IOOeTE87Tb8xsPvCFmTWJXvvgxL9XyRWOmQ00s+Vm9mszWx3FdrqZnWpmC81svZldU8F7liykpCLp0gh4gPDLuguwFbgzYf1jwGtAO+B64MKKDuTuI4FHgf9x973cfWq06jRgItA2Wj8O2AEcDBwJnAT8v2jbHwJDo+WFwPDqvhEz6wj8G7gR2Ae4EphkZvkJm50PXAy0B/aItsHMDgCeAf4M5AMFQJG7vw6si2IscSHwUAVhzAH+YmbnmlmXCrY5EegeHfM3CUnyp8DpwDeB/YENwF+qiG8syed8WMLrnAd8B2hbzSuVDkBzoCPwe+A+4HvA0UB/4Hdm1q0ax5EsoKQiaeHu69x9krtvcffNwE2ELzWiL8VjgN+7+1fu/grwVC1eZra7P+nuu4DWwKnA6OjX/GrgduDcaNsRwB3uvszd1wO31OB1vgc87e5Pu/sud38BmBu9XokH3H2hu28FJhC+nCEkm6nu/nd33x6dl6Jo3YPRsYmuHE4mJNvynA28DPwOWGpmRWZ2TJltboje+9uEhH5etPwy4Fp3X+7u2whJfHh0lVZZfBX5U3Qet1axXYntwE3uvh14HMgD/ujum939HeBd4IhqHksynJKKpIWZtTCze83so6h4agbQNiru2R9Y7+5bEnZZlrDvPVGRy+dVFI0sS3h8ANAUWBkV43wG3Eu4ciB6zcTtP6rB2zkAOLvkuNGx+wH7JWzzacLjLcBe0ePOwOIKjvsIMMzMWhKS3svuvrK8Dd19g7tf5e69gX2BIuBJM7OEzcq+v/0T4n8iIfb3gJ3RcSqLryLLqt4kyTp33xk9LklEqxLWb6X0fEmWU1KRdLkC6Akc5+6tgQHRcgNWAvuYWYuE7TuXPHD3y6Iil73c/eZKXiOx0n4ZoQI7z93bRrfW0Zcw0Wt2Tti+bBHSF0BiPImV4MuAhxOO29bdW7r7rZXElrjvQeUG774CmA2cSSj6ergax8Pd1wJjCEljn4RVZd/fJwkxnFIm/ubR61cYH8nnt7LlW6j43EkDo6QiqdDUzJon3JoArQi/QD+LinauK9nY3T8iFB9db2Z7mNkJwLByj1xN0S/854HbzKx1VDl9kJl9M9pkAvAzC01z9wauKnOIIuBcM2tqZmXrXEquKE6OKvybRxXQnajao8BgMxsRVWq3M7OChPUPAb8mtOb6Z0UHMbP/NrPDomO0Ai4HFrn7uoTNfhddIfYm1O+Mj5bfA9wU1Z9gZvlmdlo14lsFHFiN91gEnB+dmyFExZzSMCmpSCo8TUggJbfrgTuAPQl9S+YAz5bZ5wLgBEJl9Y2EL8BtdYzj+4RK8ncJldETKS2iug94jtAs9w12/wL/HeEX+wbgBhLqNtx9GaFRwDXAGsKv+19Rjf8fd/+YUPdyBbCe8AWcWH/wBFHxVJniwLJaRNt+BiyJ9vlumW1eAhYB04Ax7v58tPyPhDqr581sM+HvcVw14vsrcGhUbPZkJbH9nPCj4DPC37WybSXHmSbpkkxgZuOB9939uio3Ts3rdQWWAk3j7mthZouBSxNatdV0/65kyHsR0ZWKxCLqp3FQVEw1hHAl8GTMYdU7MzuLUEfxn7hjEUmF2HrUm1lnQnnyvoR/qrHu/seo/H080BX4EBjh7hviilPSpgOhCKodsBy43N3fjDek+mVm04FDgQujZtEiWS+24i8z2w/Yz93fiCoe5xE6aI0kNDe91cyuAvZ299/EEqSIiNRIbMVf7r7S3d+IHm8mtJ3vSCgGeTDa7EFCohERkSyQERX1UUXjDOAw4GN3bxstN2BDyfMy+4wijFxLy5Ytj+7Vq1d9hZvTitcVA9CzXc+YIxGRdJs3b95ad8+vesvqiz2pmNlehKaQN7n7P83ss8QkYmYb3H3vyo5RWFjoc+fOTXOkDcPAcQMBmD5yeqxxiEj6mdk8dy9M5TFjbf1lZk2BScCj7l7Sb2BVVN9SUu+yOq74RESkZuJs/WWEzlXvufv/Jax6ijDM+a3R/b9iCK/B+u2A38YdgohksTgn6foGYbyjt82sKFp2DSGZTDCzSwiD4o2IJ7yGafCB9TGJo4jkqtiSSjTcuVWwelB9xiKlij4tAqCgQ0GscYhIdtJ0wpJk9LOjAVXUi0jtaJgWERFJGSUVERFJGSUVERFJGSUVERFJGVXUS5KbB1U2e6+ISOWUVCRJ38594w5BRLKYir8kyaxls5i1bFbcYYhIltKViiS5Zto1gPqpiEjt6EpFRERSRklFRERSRklFRERSJjeTijt8+WXcUYiINDi5V1G/fTt07QoXXww33hh3NFnnjiF3xB2CiGSx3EsqTZtCx47w8stxR5KVNOS9iNRFbhZ/9e8Pr74K27bFHUnWmbpkKlOXTI07DBHJUrmbVLZtg7lz444k69w440ZunKFiQxGpndxMKv36hXsVgYmI1KtYk4qZ/c3MVpvZgoRl15vZCjMrim6n1vjAeXlwyCFKKiIi9SzuK5VxwJBylt/u7gXR7elaHbl/f5g5E3burEt8IiJSA7EmFXefAaxPy8H794eNG+Htt9NyeBER2V2mNin+iZl9H5gLXOHuG2p8hAEDwv3LL0NBQSpjy2n3Dr037hBEJIvFXfxVnruBg4ACYCVwW3kbmdkoM5trZnPXrFmz+wZduoSb6lVqpGdeT3rm9Yw7DBHJUhmXVNx9lbvvdPddwH3AsRVsN9bdC929MD8/v/yD9e8fkop7+gLOMZOLJzO5eHLcYYhIlsq4pGJm+yU8PQNYUNG2VerfHz79FBYvrnNcDcVts2/jttnlXhyKiFQp1joVM/s7MBDIM7PlwHXAQDMrABz4ELi01i/Qv3+4f/llOPjgOsUqIiJVizWpuPt55Sz+a8pe4JBDoF27kFQuvjhlhxURkfJlXPFXSpmF3vUzZsQdiYhIg5DbSQVCEdjixbByZdyRiIjkvEztp5I6ifUqI0bEG0sWePiMh+MOQUSyWO5fqRx5JLRsqf4q1dS5TWc6t+kcdxgikqVyP6k0bQonnKCkUk3jF4xn/ILxcYchIlkq95MKhCKw+fPhs8/ijiTj3T33bu6ee3fcYYhIlmo4ScUdZs2KOxIRkZzWMJLKcceFYjAVgYmIpFXDSCotWsDRR6u/iohImjWMpAKhCOz112Hr1rgjERHJWQ0rqWzfDq+9FnckGW3iiIlMHDEx7jBEJEs1nKTSr18YtkX1KpXKa5FHXou8uMMQkSzVcJLK3nvDYYcpqVRhXNE4xhWNizsMEclSDSepQCgCmzULduyIO5KMpaQiInXR8JLK55/DW2/FHYmISE5qeEkF6l4EtnkzvP02LFkCu3bVPS4RkRyR+6MUJ+rYEbp1C/1VRo+ufNslS6C4GJYuhQ8/DPclt/XrS7dr2RIOPzzc+vQpfbzPPrWPc8cO+M9/YPx42LQpNDIYMCAcv3Hj2h9XRCTNGlZSgXC18swzYdgWs93Xv/YaXHcdPPts6bI99oCuXUNCKiwM9127ll6xzJ8PkybBffeV7tOxIxQUhMEs+/aFY48NCagi7jBnDvz97yGZrF4NrVuHBgYToya+rVuHBNO/f0gyhYUhNoCdO8OcMR9/DMuWhfuPP4YVK0KsAwaE/dq1q9v5ExGphLl73DHUWWFhoc+dO7d6G99/P/zwh/Dee9CrV+nyN94IyWTKlPDFe8UV4Yu4Wzfo0AEaVVFS6B6+1OfPD4nm7bdh3jx4992wvnHjkGT69i29dekCCxbAY4/B44+Hq6DmzWHoUDj/fDjllPD8449Dkd3LL4errPfeC8fcc0849FBYsyYkj507k2Nq1Soktw8/hC+/DMsOPxy++c3w3gYMgH33Tdply/YtALRo2qJ651NEspaZzXP3wpQeM86kYmZ/A4YCq939sGjZPsB4oCvwITDC3TdUdpxW3Vr50dcdXelrDe0xlCv7XgkLFzLwlp6MLBjJyJ8/wNpXX2T4P86GdeugSRPo3Dl8EZcpZhpZMJKRBSNZu2UtwycM54oTrmBYz2EUry3m0imXVvzCO3bApk1c8dkhDJu5luKFs7l08FZungZ9P9+bWXtt4JpBhCuSffeFvLxyi7huHnQzfTv3ZdayWVzz7K+4d6/z6PnaYiZ/Mp3bun0CzZpDs2YhCTVrFm5NogvRXbu4t+Nl9Jy7lMnzJ3LbXvOZ+Pgu8rbAuFM6MO6EFuF9V2LiiInktcj7unXY9JHTARgzawxTFk6pdF8gafvZy2czacQkAK6eejWzl8+udN92Ldolbb9u6zrGDhsLwKjJo1i4bmGl+/do1yNp+3Z7tuOWwbcAcNaEs1i3ZV2l+5/Q6YSk7U/odEL4LAEDxw2sdF9I+OxF25f9LFWl1p+9SNntkz5L066pcv+y29879F565vVkcvFkbpt9W5X7l92+7GepKvrspe+z99LFL6U8qcRd/DUOuBN4KGHZVcA0d7/VzK6Knv8mZa/YvXsYXPK112D4cHhmEpzXOBQRdexY+kWcKk2ahPqV73wP/msYrHoH/nER9BkACzbC4W2h5ZwQU3U1bQrf/jac9xMongxV/WM3ahTGPjv5fCjuA7PGwMhrYPZbMP32cCUUJZVPNn8CwP6t9q/lGxaRhiz24i8z6wpMSbhSKQYGuvtKM9sPmO7uPSs7Ro2KvyAkk0mTQvHQL34Rbm3b1vo9ZLXrr4cbboCvvoKmTb/+9VPy605Eclc6ir/ivlIpz77uvjJ6/Cmwb3kbmdkoYBRAly5davYKv/89HH88/OAHdWullQvy88P9unWh7khEpA4yMal8zd3dzMq9lHL3scBYCFcqNTpwnz7hJqVJZc0aJRURqbNM7Py4Kir2IrpfHXM8ua19+3C/Zk28cYhITsjEpPIUcFH0+CLgXzHGkvtKrlRWK3eLSN3FWvxlZn8HBgJ5ZrYcuA64FZhgZpcAHwEj4ouwAUgs/kIV9CJSN7EmFXc/r4JVg+o1kIasXbswsoCKv0QkBTKx+EvqU+PGoQVclFTGzBrDmFljYg5KRLKVkoqEIrAoqUxZOKVavZRFRMqjpCKhBZiKv0QkBZRUJOlKRUSkLpRUJCQVNSkWkRTI6B71Uk/y88PEYzt3smfTPeOORkSymJKKhKTiDuvW8cwFz8QdjYhkMRV/yW4dIEVEaktJRZKSyh9e+gN/eOkP8cYjIllLSUWSBpWctnQa05ZOizceEclaSiqi4i8RSRklFQnjf4GaFYtInSmpSJjzfu+9daUiInWmJsUSRL3q27VoF3ckIpLFlFQkiJLKpBEvxh2JiGQxFX9JoPG/RCQFlFQkiEYqvnrq1Vw99eq4oxGRLKXiLwny82HtWmYvnx13JCKSxTI2qZjZh8BmYCeww90L440ox+Xnw65dsH17aA0mIlILGZtUIie6+9q4g2gQSjpAKqmISB2oTkWCxKQiIlJLmZxUHHjezOaZ2aiyK81slJnNNbO5a9Rqqe6ipNLJW9GpdaeYgxGRbJXJxV/93H2FmbUHXjCz9919RslKdx8LjAUoLCz0uILMGdGgko80Px/OvCzmYEQkW2XslYq7r4juVwNPAMfGG1GOy8sL97rqE5E6yMikYmYtzaxVyWPgJGBBvFHluD32gDZtGP3FJEY/OzruaEQkS2Vq8de+wBNmBiHGx9z92XhDagDy8yna9Ql8WhR3JCKSpTIyqbj7EuCIuONocPLzYbtacItI7WVk8ZfEJD9fTYpFpE4y8kpFYpKfD9u/ijsKaQjWr4f586FNGzjyyLijkRRSUpFS7dvT450d0K573JFIrti+HYqLQwJJvK1YEda3bg3r1kETfRXlCv0lpVR+PmOfchj3P3FHItnms89C8li4MPn+/ffhq+jqt2lTOOQQOPFEOOII2LABbr4ZioqgUEP75QolFSlVMlTLmjVhemGRstxh0SJ46SWYMyckjuLi5P5NjRvDgQdC9+4wZAj06RNuPXuGpuslli8PSeWVV5RUcoiSipTKz2fUMOClKxnb46m4o5FM4B6uNl56qfS2cmVY164d9O4Np50WEkaPHuG+W7fk5FGRTp2ga9eQVEaPTue7kHqkpCKl8vNZ2A7YuCTuSCQO7iFhvP12uL36KsyYAatXh/X77Qff/GbprVcvCH3Jaq9/f3j++fDadT2WZAQlFSmlkYobjs2bYcGC0gRSclu/vnSbLl3g5JNLk8hBB6X+i79fP3j44VCk1l0NRHKBkoqUKkkqX6lZcc7Yvj1UmpdNHh9+WLpNy5Zw2GFw1llw+OHhdthhpePBpVO/fuH+lVeUVHKEkoqUat48VLLqSiV7rV0LU6bAtGmh6W5i66vGjUO9x7HHwiWXlCaQrl2hUUz9oHv1gn32CUnl4ovjiUFSSklFkhRsbgl77RV3GFITixfDv/4Vbq+8EqaF3nff0Knw5JNLk8chh0CzZnFHm6xRo3C18sorcUciKaKkIknu+OgQ2NA67jCkMu4wbx48+WRIJAuiAbwPPxyuvTa0xjrqqOyp+O7XD556KjQIiOb1keylpCLJ8vNh2bK4o5BEO3eGepAZM+Dll0tbZDVqFFpP/d//hURy4IFxR1o7ifUqZ54ZbyxSZ0oqkuR7vd6FNit5JO5AGrKvvgpXIjNmhNvMmbBxY1h3wAFw0kkweDB85zv1U5mebkcfHerzlFRygpKKJFm+53bY40v1G6hvH38MTz8dbtOmwZYtYXmvXnDOOeGKpH//kFRyzR57wHHHqV4lR1QrqZjZPuUs3uzuaiaUa5o2DQll06Ywgqykx44dMGtWSCL//ndpvUi3bqEV1Le+FYqFGkodQ79+cOut8MUXoYmzZK3qXqm8AXQGNgAGtAU+NbNVwA/dfV56wpN6VzK8xpo1SiqptmlTSCBPPgnPPReKtJo0gQEDYMyYUJzVs2fDvELs1y/UHc2ZA4MGxR2N1EF1k8oLwER3fw7AzE4CzgIeAO4CjktPeFLvmjYN92vWwMEHxxtLLvjsM5g8GSZODIlk2zbo0CF0NDz1VPj2t8Pw7w3dCSeEZPrKK0oqWa66SeV4d/9hyRN3f97Mxrj7pWaWlobvZjYE+CPQGLjf3W9Nx+tIshPaHwXL3k4edVZqZv360ER24sQwrtX27WHwxMsvh+HDwxdoXJ0NM1WbNmE4fNWrpMeXX8KSJaEzbOK8NmlQ3aSy0sx+AzwePT8HWGVmjYFdqQ4qOu5fgG8Dy4HXzewpd3831a8lyW7pfwNc+CCcuzruULLLtm2hz8i4cfDCC6HO5IAD4Gc/C4nk2GOVSKrSrx888EA4d5q0q+a+/DJ0hP3ggzCWWsn9okWhm4B72K5RozAkzlFHhe1TrLp/ufOB64Ano+czo2WNgREpjwqOBRa5+xIAM3scOA1QUkm3xDlVpHLu8MYb4YvwscfCpFOdOsEvfwlnnx2ayjbE+pHa6tcP7rxTk3ZVh3sYv2327FAPNWdOOG+JQyy1axeSx4AB4f7gg8MwPYceCi1ahG3S8PmsVlJx97XATytYvSh14XytI5DYA285ZeptzGwUMAqgS5cuaQihYTpryoVwXmMmKalUbM0aePRR+NvfQqfEZs3gjDNCq61Bg8IYW1Jz3/hGuNekXbvbsCEkjVdfLU0kJVMStGgBxxwDV1wRJkPr3j2MKB3TRHvVbVKcD/wa6A00L1nu7t9KU1xVcvexwFiAwsJCjyuOXLNuyzpo00RXKmXt2hXqR8aODRXvO3aEf+S77oJzz9VMmamgSbvC5+qDD+Ctt5LrPhJHuejRI8yoefzxoX7usMMyqriwupE8CowHhgKXARcB6fzWWUFowlyiU7RM6kPTPZRUSqxaFYq3xo6FpUtD8eDPfw4jR4Z/ZkmthjRp15o1yYlj/nx4551QPwchUfTqFc5Jnz6hIcMxx4RirQxW3aTSzt3/amY/d/eXgJfM7PU0xvU60N3MuhGSybmEOhypD02bNuyk4h6mzb3nHvjnP0M59cCBcMstoZirOlPlSu3k4qRdX321e6ur+fNLp2WG0My8Tx/46U/DfZ8+IaFk2qjS1VDdpFJS+7PSzL4DfAKU18s+Jdx9h5n9BHiO0Bjgb+7+TrpeT8poqEll/Xp46KGQTIqLoW1b+PGP4dJLwz+4pF82T9rlDp9+Wpo0Soqw3nsvFGtB+EHSu3cYv60kefTpk1MjJ1Q3qdxoZm2AK4A/A62B0ekKCsDdnwaeTudryO4GdRsEC3fA6rkNowjCPVR+3nMPjB8fmmUef3xoGjxiBOy5Z9wRNizZNGmXO7z7bqhjmzo1JJG1a0vXd+oUEsZ3vlNafNW9e2kH4xxV3dZfU6KHG4ETAcxsdJpikhj97pu/g1ebwbaZ8Pnn0KpV3CGlx+bNoQXXPfeEL4O99gr1JJdeCgUFcUfXcGX6pF1ffRVGjp48OdyWLg3LCwrg9NNLrzwOPzwkxwaoLk0GfgnckaI4JJMk9lXJtaRSVBQSyaOPhqR5xBHh+fnn5957zVYlk3atWhVmsIzbypXhSmTy5DDUzqZNYaj+QYPgqqvClUjHjnFHmTHqklRyvFykYTrl0VPgszU8AyGpZOvET4m2b4dJk+CPfwzt+5s3D82AL7ss9HTP9SK+bFNSrzJzZv3Or+IerjzefDN0ai25X7UqrO/QIRSJDhsW5rMp6UAoSeqSVNQ3JAdt3b4VmkR/2myvrN+wITQFvvNOWL48lGfffjtcdJH6lWSy+pq0a+vWUJQ1dSrMnRuSSMlkaI0bhwr1IUPgyCOhb98Ql4baqVKlScXMNlN+8jBANZi5KnGk4mxUXAx/+lOobN+yJcxNcvfdYVRgfSlkvnRN2uUePhvPPhuKsaZPDw0zmjULdSLnnRcSyFFHhT5IzZtXdUQpR6VJxd1VyNwQZWNScQ8zJt5xR5izZI894IILQs/sPn3ijk5qqmTSrs8/D40oamvz5nAlUpJIPvooLO/ZMzTKGDIkjI2loqyUyZy+/ZI5GjcOv9KyJals2gSXXBKGmm/fHq6/PtSXZEIlr9ROyaRdr75a8/lVPvqotHXWiy+GOrVWrcJxrr4aTj45DAcjaaGkIkmG9hgaHuQvKR2wLpO9/XaY8GrJktDjffRoFVvkgppM2rVrF7z+emkiKZknpEePMKTO0KGhTiTH+4dkCiUVSXJl3yvDg/aPZ/6VyoMPhomv2rYNv0j79487IkmVyibt2rUr9FKfOTPcnnsutNBq3Dhc4YwZE1po9ehR/3GLkopUID8/c5PK1q1h8qv774cTTwxzmXToEHdUkmolk3Zt3Ajz5sGsWeE2e3aYphnC5/TEE+G734VTTmmwHQ4ziZKKJBk4biAA0/O7hF+DmWbx4jCTYlERXHMN3HBDRg37LSlUMmnX3nuXzlrYu3eYAO0b3whFWgcfrH5GGUb/jVK+TLxSefLJMJRKo0YwZUroySy565RT4MILoUuXkESOP179i7KAkoqULz8/9PHYsiX+5paffAI33hj6mhQWwj/+odY7DUHr1mHUaMkq6gkm5cuEuepXrAh1JwceGHrG/+xnoeJWCUUkYympSPlKkkoczYpLkslBB4Xpei+4ABYuDGN3ZeGkRSINiYq/JMmI3iPCg53RpEH1eaWyYkXoRX3ffaHj20UXhcr4XBjUUqSBUFKRJD865kfhweLF4b4+ksqWLaGn8z33hD4II0eGZNKtW/pfW0RSSklFkmzZvgWAFvVVp7JiRehj8OabYaiVa69VnYlIFsu4pGJm1wM/BEq+za6JphaWenDqo6cCMP2iF8OgjOlMKvPmhYSyaVOYlGno0PS9lojUi4xLKpHb3X1M3EE0aGbp7asycSJ8//vhNWbO1EjCIjlCrb+kYulIKu5w002hV3RBAbz2mhKKSA7J1KTyEzObb2Z/M7Nyu9Ca2Sgzm2tmc9dkWs/vXJGfn9omxV9+Ga5Ofvvb0Ez4P//R8PQiOSaWpGJmU81sQTm304C7gYOAAmAlcFt5x3D3se5e6O6F+SWVypJa7dun7kpl9eowhPkjj4Te8Q8/rCHqRXJQLHUq7j64OtuZ2X3AlDSHIwlGFowsfZKq4q933w1T+a5eHYZYGT687scUkYyUcRX1Zrafu6+Mnp4BLIgznoZmt6Ty+eeh2Kq2VxXvvx+GJm/UCGbMCGN3iUjOyrikAvyPmRUADnwIXBprNA3M2i1rAchrkZc8/lfnzjU/2AcfwLe+FVqSvfgi9OqVwkhFJBNlXFJx9wvjjqEhGz4hFE1NHzm9bkllyZKQULZvh+nTlVBEGoiMSyqSQWrbq/6jj0JC2bIltPDq3Tv1sYlIRlJSkYq1jwaVrEmz4uXLQ0LZuBGmTQvzjItIg6GkIhWr6ZXKypUhoaxZA1OnwlFHpS82EclISipSsTZtoGnT6iWVVatCQvnkE3j+eTj22PTHJyIZR0lFklxeeHnpEzPIy6s6qaxdC4MHh7qUZ5+Fvn3TG6SIZCwlFUlyzmHnJC+oqgPk4sVw+umwaBFMmQIDBqQ1PhHJbJk69pfEZNnGZSzbuKx0QWVJ5dlnQ2fGFStCQhk0qH6CFJGMpaQiSS584kIufCKhq1B5SWXXrjDS8KmnQpcuMHeuEoqIACr+kqq0b5/cpHjTpjB3/JNPwvnnw9ix0LJlbOGJSGZRUpHK5eeHRLJtGyxdCmecEYZfuf12+PnPQ2W+iEhExV9SuZK+KvffH5oJr1sX+qCMHq2EIiK7UVKRypUklZ/8BA45BN54AwYOjDUkEclcKv6SJFeccEXygp49wxXJJZfAn/+sibVEpFJKKpJkWM9hyQt694b166Ft21jiEZHsouIvSVK8tpjitcXJC5VQRKSadKUiSS6dEuZEmz5yeryBiEhW0pWKiIikjJKKiIikTCxJxczONrN3zGyXmRWWWXe1mS0ys2IzOzmO+EREpHbiqlNZAJwJ3Ju40MwOBc4FegP7A1PNrIe776z/EEVEpKZiSSru/h6A7d4j+zTgcXffBiw1s0XAscDs+o2w4frtgN/GHYKIZLFMa/3VEZiT8Hx5tGw3ZjYKGAXQpUuX9EfWQAw+cHDcIYhIFktbUjGzqUCHclZd6+7/quvx3X0sMBagsLDQ63o8CYo+LQKgoENBrHGISHZKW1Jx99r85F0BdE543ilaJvVk9LOjAfVTEZHaybQmxU8B55pZMzPrBnQHXos5JhERqaa4mhSfYWbLgROAf5vZcwDu/g4wAXgXeBb4sVp+iYhkj7hafz0BPFHBupuAm+o3IhERSYVMK/4SEZEslmlNiiVmNw+6Oe4QRCSLKalIkr6d+8YdgohkMRV/SZJZy2Yxa9msuMMQkSylKxVJcs20awD1UxGR2tGVioiIpIySioiIpIySioiIpIySioiIpIwq6iXJHUPuiDsEEcliSiqSREPei0hdqPhLkkxdMpWpS6bGHYaIZCldqUiSG2fcCGgGSBGpHV2piIhIyiipiIhIyiipiIhIyiipiIhIyqiiXpLcO/TeuEMQkSwW1xz1Z5vZO2a2y8wKE5Z3NbOtZlYU3e6JI76GrGdeT3rm9Yw7DBHJUnFdqSwAzgTK+1m82N0L6jccKTG5eDIAw3oOizkSEclGsSQVd38PwMzieHmpxG2zbwOUVESkdjKxor6bmb1pZi+ZWf+4gxERkepL25WKmU0FOpSz6lp3/1cFu60Eurj7OjM7GnjSzHq7+6Zyjj8KGAXQpUuXVIUtIiJ1kLak4u41HufD3bcB26LH88xsMdADmFvOtmOBsQCFhYVet2hFRCQVMqr4y8zyzaxx9PhAoDuwJN6oRESkumKpqDezM4A/A/nAv82syN1PBgYA/2Vm24FdwGXuvj6OGBuqh894OO4QRCSLxdX66wngiXKWTwIm1X9EUqJzm85xhyAiWSyjir8kfuMXjGf8gvFxhyEiWUrDtEiSu+feDcA5h50TcyQiko10pSIiIimjpCIiIimjpCIiIimjpCIiIimjinpJMnHExLhDEJEspqQiSfJa5MUdgohkMRV/SZJxReMYVzQu7jBEJEspqUgSJRURqQslFRERSRklFRERSRklFRERSRklFRERSRk1KZYkT1/wdNwhiEgWU1KRJC2atog7BBHJYir+kiR3vX4Xd71+V9xhiEiWUlKRJBPemcCEdybEHYaIZCklFRERSZlYkoqZ/a+ZvW9m883sCTNrm7DuajNbZGbFZnZyHPGJiEjtxHWl8gJwmLv3ARYCVwOY2aHAuUBvYAhwl5k1jilGERGpoViSirs/7+47oqdzgE7R49OAx919m7svBRYBx8YRo4iI1FwmNCn+ATA+etyRkGRKLI+W7cbMRgGjoqfbzGxB2iJMnTxgbdxBVEOeXWxZESdZcj5RnKmSDTFC9sTZM9UHTFtSMbOpQIdyVl3r7v+KtrkW2AE8WtPju/tYYGx0nLnuXliHcOuF4kwtxZla2RBnNsQI2RVnqo+ZtqTi7oMrW29mI4GhwCB392jxCqBzwmadomUiIpIF4mr9NQT4NfBdd9+SsOop4Fwza2Zm3YDuwGtxxCgiIjUXV53KnUAz4AUzA5jj7pe5+ztmNgF4l1As9mN331mN441NX6gppThTS3GmVjbEmQ0xQgOO00pLnkREROpGPepFRCRllFRERCRlMjKpmNmQaJiWRWZ2VTnrm5nZ+Gj9q2bWNWFducO8VHXM+ozTzL5tZvPM7O3o/lsJ+0yPjlkU3drHFGNXM9uaEMc9CfscHcW+yMz+ZFHFWExxXpAQY5GZ7TKzgmhdSs9lNeMcYGZvmNkOMxteZt1FZvZBdLsoYXkc57PcOM2swMxmm9k7FoZROidh3TgzW5pwPgviijNatzMhlqcSlneLPiOLos/MHnHFaWYnlvl8fmlmp0fr4jifvzSzd6O/7TQzOyBhXWo+n+6eUTegMbAYOBDYA3gLOLTMNj8C7okenwuMjx4fGm3fDOgWHadxdY5Zz3EeCewfPT4MWJGwz3SgMAPOZVdgQQXHfQ04HjDgGeCUuOIss83hwOJ0nMsaxNkV6AM8BAxPWL4PsCS63zt6vHeM57OiOHsA3aPH+wMrgbbR83GJ28Z5PqN1n1dw3AnAudHje4DL44yzzGdgPdAixvN5YsLrX07p/3vKPp+ZeKVyLLDI3Ze4+1fA44ThWxKdBjwYPZ4IDIqyZ0XDvFTnmPUWp7u/6e6fRMvfAfY0s2Z1jCelMVZ0QDPbD2jt7nM8fOIeAk7PkDjPi/ZNlyrjdPcP3X0+sKvMvicDL7j7enffQBj/bkhc57OiON19obt/ED3+BFgN5NcxnpTHWZHoM/EtwmcEwmfm9AyJczjwjCd3o0il6sT5YsLrJw6RlbLPZyYmlY7AsoTn5Q3V8vU2HsYQ2wi0q2Tf6hyzPuNMdBbwhrtvS1j2QHQ5/Ls6FoXUNcZuZvammb1kZv0Ttl9exTHrO84S5wB/L7MsVeeyunHWdN+4zmeVzOxYwi/exQmLb4qKTm5PwQ+husbZ3MzmmtmckiIlwmfiMy8dWzBjzifhCrvs5zPO83kJ4cqjsn1r/PnMxKTSYJhZb+C/gUsTFl/g7ocD/aPbhXHERij26OLuRwK/BB4zs9YxxVIlMzsO2OLuiWPAZcq5zDrRL9SHgYvdveTX99VAL+AYQjHJb2IKr8QBHoZCOR+4w8wOijmeCkXn83DguYTFsZ1PM/seUAj8b6qPnYlJpTpDtXy9jZk1AdoA6yrZNx3Dv9QlTsysE/AE8H13//qXoLuviO43A49Rt1Gaax1jVIS4LoplHuHXao9o+04J+8d+LiO7/QpM8bmsbpw13Teu81mh6MfDvwnj9H09wKu7r/RgG/AA8Z7PxL/vEkL92ZGEz0Tb6DNS42OmI87ICOAJd99esiCu82lmg4FrCSOabKti35p/PlNVSZSqG6GX/xJCRXtJZVPvMtv8mORK2wnR494kV9QvIVReVXnMeo6zbbT9meUcMy963JRQLnxZTDHmA42jxwdGH6R9vPyKu1PjOpfR80ZRfAem61xWN86Ebcexe0X9UkIl6N7R49jOZyVx7gFMA0aXs+1+0b0BdwC3xhjn3kCz6HEe8AFRpTTwD5Ir6n8UV5wJy+cAJ8Z9PgmJdzFRY4x0fD5r/QbSeQNOJUzetZjwawngvwiZFaB59MFZFL3hxC+Ta6P9iklopVDeMeOKE/gt8AVQlHBrD7QE5gHzCRX4fyT6Yo8hxrOiGIqAN4BhCccsBBZEx7yTaGSGGP/mAwlD/SQeL+XnsppxHkMod/6C8Kv5nYR9fxDFv4hQrBTn+Sw3TuB7wPYyn82CaN1/gLejWB8B9ooxzr5RLG9F95ckHPPA6DOyKPrMNIv5796V8KOnUZljxnE+pwKrEv62T6X686lhWkREJGUysU5FRESylJKKiIikjJKKiIikjJKKiIikjJKKiIikjJKK5JQyI9cWWcII1tnOzI40s79Gj0ea2Z1l1k83s8JK9n/czLqnO05p2OKaTlgkXba6e0F5K6Kxv8xLhx3JNtcAN9Zh/7uBXwM/TE04IrvTlYrkNAvzwhSb2UOEDlydzexXZvZ6NJDfDQnbXmtmC83sFTP7u5ldGS3/+grAzPLM7MPocWMz+9+EY10aLR8Y7TPRzN43s0dLBrM0s2PMbJaZvWVmr5lZKzObkTiXRvT6R5R5H62APu7+VjXe83cTrtSKzWxptOplYHDCECYiKacPl+SaPc2sKHq8FPgF0B24yN3nmNlJ0fNjCcNOPGVmAwg9oc8FCgj/F28QeuRX5hJgo7sfE40wO9PMno/WHUkYNugTYCbwDTN7DRgPnOPur0djbG0F/gqMBEabWQ+geTnJo6RXc6JzzKxfwvODAdz9KeApADObALwULd9lZouAI6rx3kRqRUlFck1S8VdUp/KRlw6MeFJ0ezN6vhchybQiDPi3Jdrv65kEK3ES0MdKZ/prEx3rK+A1d18eHauIMFTHRmClu78O4O6bovX/AH5nZr8iDJUxrpzX2g9YU2bZeHf/ScJ7nZ640sx+TTgff0lYvJow+ZaSiqSFkoo0BF8kPDbgFne/N3EDMxtdyf47KC0qbl7mWD9198ThzDGzgUDi/Dg7qeR/zd23mNkLhAmVRgBHl7PZ1jKvXaloJNqzgQFlVjWPjiWSFqpTkYbmOeAHZrYXgJl1tDB3/QzgdDPbM6q/GJawz4eUftEPL3Osy82saXSsHmbWspLXLgb2M7Njou1bJdRv3A/8CXjdw8x7Zb1HVLxVFQvzjv8FONvdyyaQHuxejCaSMrpSkQbF3Z83s0OA2VHd+efA99z9DTMbTxj1djXwesJuY4AJZjaKMM9IifsJxVpvRBXxa6hkqlV3/8rMzgH+bGZ7Eq4YBhPmWp9nZpsI82qUt+/7ZtbGzFp5mB+mMiMJMyA+Gb3HT9z9VDPbl1Ac9mkV+4vUmkYpFimHmV1P+LIfU0+vtz9hoqleFTV5NrNfAJvd/f5avsYvgE3u/tdaBypSBRV/icTMzL4PvEqY/6KyPjR3k1xXU1OfAQ/WYX+RKulKRUREUkZXKiIikjJKKiIikjJKKiIikjJKKiIikjJKKiIikjL/H3DbH9bCFmNTAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.figure()\n", "\n", "# Plot lag-frequency spectrum.\n", "plt.plot(cross.freq, lag, 'r')\n", "\n", "# Find cutoff points\n", "v_cutoff = 1.0/(2*delay)\n", "h_cutoff = lag[int((v_cutoff-0.0050)*1/0.0050)]\n", "\n", "plt.axvline(v_cutoff, color='g',linestyle='--')\n", "plt.axhline(h_cutoff, color='g', linestyle='-.')\n", "\n", "# Define axis\n", "plt.axis([0,0.2,-20,20])\n", "plt.xlabel('Frequency (Hz)')\n", "plt.ylabel('Lag')\n", "plt.title('Lag-frequency Spectrum')\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "According to Uttley et al (2014), the lag-frequency spectrum shows a constant delay until the frequency (1/2*time_delay) which is represented by the green vertical line in the above figure. After this point, the phase wraps and the lag becomes negative. " ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "13it [00:00, 3156.72it/s]\n" ] } ], "source": [ "cross = AveragedCrossspectrum(lc, var, segment_size=200)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZUAAAEWCAYAAACufwpNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAvIUlEQVR4nO3debyUdfn/8dcFKJssyiFF2bSEfm5hHO2rolL6U/KruYYaWqSFWdaXzA2wb1ailpqUqUmalJBCmiZqZtgPrYAU7GhuEG4skiKioiDbuX5/fO7juWeYOWfOOfece5b38/GYx5m512tuhrnms97m7oiIiCShQ9oBiIhI5VBSERGRxCipiIhIYpRUREQkMUoqIiKSGCUVERFJjJKKlAUz29nMHjOzdWZ2bdrxiEhuSirSJmb2ipkd2Q6nGge8CfR09++0w/lKipn1N7O7zexNM3vHzJ4xs7FFPuc0M7u8mOeQytMp7QBECjQIeM7zjNY1s07uvqWdY2pPtwNPEa7DRmBfYJc0A6qCay6toJKKFIWZ7Whm95vZajNbGz3vH1u/e6w6a46Z3WBm0/McaxrwJeAiM3vPzI40s8vM7C4zm25m7wJjzayXmd1qZqvMbKWZXW5mHaNjdDSza6Jf+i+Z2TfMzM2sU7Q+o8QVHX967PV/mdk8M3vbzJ4ys5GxdXPN7Idm9vfo/TxsZjWx9SNi+y43s7FmdoCZvd4QX7TdSWb2VJ5LegAwzd3fd/ct7v5Pd/9jtN/g6L2MM7PXovd/Qey4HczsEjN70czWmNksM9upmfjGAWNi13x27DpdbGZPA++bWafo3B+L/3s1lHDMbKSZrTCzi8zsjSi2E8zsGDNbYmZvmdnEPO9ZypCSihRLB+A2wi/rgcAG4Oex9b8FHgf6AJcBZ+Y7kLuPBWYAP3b3Hdx9TrTqeOAuoHe0fhqwBfgYsD9wFPCVaNuvAsdGy2uBUwp9I2a2G/AAcDmwE3ABcLeZ9Y1t9gXgy8BHgO2jbTCzQcAfgeuBvsAwoM7dnwDWRDE2OBP4TZ4wFgA3mNlpZjYwzzafBvaMjnlxLEl+EzgBOBzYFVgL3NBMfFPJvObHxc5zOvDfQO8CSyq7AF2A3YD/BX4JnAEMBw4FvmtmuxdwHCkDSipSFO6+xt3vdvf17r4OmEz4UiP6UjwA+F933+TufwPua8Vp5rv7ve5eD/QEjgHGR7/m3wCuA06Lth0NTHH35e7+FnBlC85zBvCguz/o7vXu/mdgYXS+Bre5+xJ33wDMInw5Q0g2c9z9DnffHF2Xumjdr6NjE5UcjiYk21w+D/wV+C7wspnVmdkBWdt8P3rv/yIk9NOj5V8DJrn7CnffSEjip0SltKbiy+dn0XXc0Mx2DTYDk919M3AnUAP81N3XufuzwHPAJwo8lpQ4JRUpCjPrZmY3m9mrUfXUY0DvqLpnV+Atd18f22V5bN9fRFUu7zVTNbI89nwQsB2wKqrGeRu4mVByIDpnfPtXW/B2BgGfbzhudOwRQL/YNv+JPV8P7BA9HwC8mOe404HjzKw7Ien91d1X5drQ3de6+yXuvjewM1AH3GtmFtss+/3tGov/nljszwNbo+M0FV8+y5vfJMMad98aPW9IRK/H1m+g8XpJmVNSkWL5DjAU+JS79wQOi5YbsArYycy6xbYf0PDE3b8WVbns4O5XNHGOeKP9ckIDdo27944ePaMvYaJzDohtn12F9D4QjyfeCL4cuD123N7u3t3dr2oitvi+H80ZvPtKYD5wEqHq6/YCjoe7vwlcQ0gaO8VWZb+/12IxfDYr/i7R+fPGR+b1bWr5evJfO6kySiqShO3MrEvs0QnoQfgF+nZUtfO9ho3d/VVC9dFlZra9mR0EHJfzyAWKfuE/DFxrZj2jxumPmtnh0SazgG9Z6Jq7I3BJ1iHqgNPMbDszy25zaShRHB01+HeJGqD707wZwJFmNjpq1O5jZsNi638DXETozfX7fAcxsx+Z2T7RMXoA5wJL3X1NbLPvRiXEvQntOzOj5b8AJkftJ5hZXzM7voD4Xgf2KOA91gFfiK7NKKJqTqlOSiqShAcJCaThcRkwBehKGFuyAHgoa58xwEGExurLCV+AG9sYxxcJjeTPERqj76KxiuqXwJ8I3XKfZNsv8O8SfrGvBb5PrG3D3ZcTOgVMBFYTft1fSAH/f9x9GaHt5TvAW4Qv4Hj7wT1E1VNZ1YHZukXbvg28FO3zuaxtHgWWAo8A17j7w9HynxLarB42s3WEf49PFRDfrcBeUbXZvU3E9j+EHwVvE/5dm9pWKpzpJl1SCsxsJvCCu3+v2Y2TOd9g4GVgu7THWpjZi8A5sV5tLd1/MCXyXkRUUpFUROM0PhpVU40ilATuTTmsdmdmJxPaKP6SdiwiSUhtRL2ZDSDUJ+9M+E811d1/GtW/zwQGA68Ao919bVpxStHsQqiC6gOsAM5193+mG1L7MrO5wF7AmVG3aJGyl1r1l5n1A/q5+5NRw+MiwgCtsYTupleZ2SXAju5+cSpBiohIi6RW/eXuq9z9yej5OkLf+d0I1SC/jjb7NSHRiIhIGSiJhvqoofExYB9gmbv3jpYbsLbhddY+4wgz19K9e/fhH//4x9sr3Iq2eM1iAIb2GZpyJCJSbIsWLXrT3fs2v2XhUk8qZrYDoSvkZHf/vZm9HU8iZrbW3Xds6hi1tbW+cOHCIkdaHUZOGwnA3LFzU41DRIrPzBa5e22Sx0y195eZbQfcDcxw94ZxA69H7S0N7S5vpBWfiIi0TJq9v4wwuOp5d/9JbNV9hGnOr4r+/iGF8KrWpYddmnYIIlLG0rxJ1yGE+Y7+ZWZ10bKJhGQyy8zOJkyKNzqd8KrTkXu0x00cRaRSpZZUounOLc/qI9ozFmlU9586AIbtMizVOESkPOl2wpJh/EPjATXUi0jraJoWERFJjJKKiIgkRklFREQSo6QiIiKJUUO9ZLjiiKbu3isi0jQlFclw8ICD0w5BRMqYqr8kw7zl85i3fF7aYYhImVJJRTJMfGQioHEqItI6KqmIiEhilFRERCQxSioiIpIYJRUREUmMGuolw5RRU9IOQUTKmJKKZNCU9yLSFqr+kgxzXprDnJfmpB2GiJQplVQkw+WPXQ7oDpAi0joqqYiISGJSTSpm9isze8PMnoktu8zMVppZXfQ4Js0YRUSkcGmXVKYBo3Isv87dh0WPB9s5JhERaaVUk4q7Pwa8lWYMIiKSnFJtqD/PzL4ILAS+4+5r0w6oWtx87M1phyAiZSzt6q9cbgI+CgwDVgHX5trIzMaZ2UIzW7h69ep2DK+yDa0ZytCaoWmHISJlquSSiru/7u5b3b0e+CVwYJ7tprp7rbvX9u3bt32DrGCzF89m9uLZaYchImWq5Kq/zKyfu6+KXp4IPNPU9pKsa+eHguFxQ49LORIRKUepJhUzuwMYCdSY2Qrge8BIMxsGOPAKcE5a8YmISMukmlTc/fQci29t90BERCQRJdemIiIi5UtJRUREElNyDfWSrttPvD3tEESkjCmpSIYBvQakHYKIlLHqqP4aOTI8pFkzn5nJzGdmph2GiJQplVQkw00LbwLg1H1OTTkSESlH1VFSERGRdlH5SWXGDFiwAB59FAYPDq9FRKQoKjupzJgB48bBxo3h9auvhtdKLCIiRVHZSWXSJFi/PnPZ+vVhuYiIJK6yG+qXLWvZcuGu0XelHYKIlLHKLqkMHNiy5UJNtxpqutWkHYaIlKnKTiqTJ0O3bpnLunULyyWnaXXTmFY3Le0wRKRMVXb115gx4e/ZZzc21p9+euNy2UZDQhk7bGyqcYhIearspAIhgfzyl+AOW7fC7Nmwbh306JF2ZCIiFaeyq7/izOAnP4E33oCrrip8vxkzwviWDh2KO86lvc4jIlJE1ZFU5s4NjwMPDCWXa68NY1aa0zDO5dVXQ0mnWONc2vM8SlwiUkTVkVTirrgilFomTGh+20LHuTT3Zd3c+gkT2j6eppAY2iNxiUhVM3dPO4Y267F7Dx/+veFNbnPskGO54OALAPjN5wbxxdnLYP583tzvY5wy65Sc+/zly4/mzLoOPH/hWPY6/Vu8tmAOO42/hC6b6j9c/8H2Hbh67BAeOWhnjpj/OhdOW5KxfnNHY33tJ+i1fQ82Lnmezq+/mfP8Dnz6V4eBGVPfHcmQn0zDly3j9Z22p/7yy9n1axcwe/Fs6n5y0Tbn2NjJuOuo3Xhh95503biVi+/6Dx3ffmebc7zXrw/HXrHPh6/rPRyjg2W+87tG30VNt5oPe4fNHTsXgGvmXcP9S+7PGX9cfPv5K+Zz9+i7AZgwZwLzV8xvct8+3fpkbL9mwxqmHjcVgHGzx7FkzZIm9x/SZ0jG9n269uHKI68E4ORZJ7Nm/Zom9z+o/0EZ2x/U/6APP0sjp41scl/I/OyNnDaSscPGMnbYWN5c/2bez15c9vbfOeg7HDf0OBa/uZhz7j+n2f2zt7/iiCs4eMDBzFs+j4mPTGx2/+ztbz72ZobWDGX24tlcO//aZvfP3j77s9QcffaK99l79MuPLnL32mYP1AKpNtSb2a+AY4E33H2faNlOwExgMPAKMNrd1yZ53juOGcgp89+h2/nnw5//kHe7d3foRO/3tmyzfEtHY6+rp8HV0+jXoQNWX5+xvsumei6ctoQTH3mNj7/8Lp0yV7PdVqfnE0/DIYfwzmEH0vX+P9Fjw9ZtzmPAzAv+wYqdu/CxF+fBps0YsMuajdSPnwT/eJa969dw1B1L6Lw58ySdtzhjHlzR7LXo/p/MD3V2MhERaRF3T+0BHAZ8EngmtuzHwCXR80uAHzV3nOHDh3uL3XKLO7jPnLntuq1b3S+9NKzv0CH8bXh06+Y+fbr7ypXuv/lN5rrsxxFH5F9n1ni+6dPDcePru3Z1/8pX3E88MWyb7zjZ+2Wfo67OfelS9/79c28zaFDGW7/h8Rv8hsdvaPn1FJGyAyz0pL/Xkz5giwMIJZJ4UlkM9Iue9wMWN3eMViWVLVvc99vPffBg9w0bGpe/+6778ceHS3P22e7TpoUvXrPwd/r0zOMMGtT0l3Vz6xtMn57/PPmSipl7fX1h58iVuBoSZMzhtx3uh992eMuvp4iUnWpJKm/Hnlv8ddZ+44CFwMKBAwe27orOmRMuQceO4e+uu7rvtlt4/bOfhS/s5jT3ZV3gl3mTmksahZ5j+nT3Tp0a980Rg5KKSPWouqQSvV7b3DFaVVJxD1+qDQkl/pgwoeXHaao009z6Qo7fXNIo9ByHHOL+mc/kPZWSikj1KEZSSb33l5kNBu73xob6xcBId19lZv2Aue4+tKlj1NbW+sKFC1t+8sGDc49XGTQIXnml5ccrphkzQhfjZcvChJiTJ7duupkTToCXX4annsq5uqFHSUOPGRGpXGaWeO+vUuzqcx/wpej5l4D83bPaqpymxh8zJiS6+vrwt7Xzl9XUwJu5uzCLiLRV2l2K7wBGAjVmtgL4HnAVMMvMzgZeBUYXLYCBA3OXVCp5avyGpOIeBoFmUQlFRNoi1aTi7qfnWXVEuwQweXIYVR4fzV7pU+PX1MCmTfDee5pUU0QSV4rVX+1nzBiYOjW0oZiFv1OnVvbU+DXRDbjyVIFdM+8arpl3TTsGJCKVpLqTCiTXVlEumkkq9y+5v6CpL0REclFSqTbNJBURkbZQUqk2SioiUkRKKtVGSUVEiqjybycsmXr1gk6dYPXqnKu7bte1nQMSkUqipFJtzJocAPnHMX9s54BEpJKo+qsaaVS9iBSJkko1aiKp/PDRH/LDR3/YzgGJSKVQUqlGTSSVR15+hEdefqSdAxKRSqGkUo1U/SUiRaKkUo1qamDNmjCLgIhIgpRUqlFNTUgob7+ddiQiUmHUpbgaxQdA7rRTxqo+3fqkEJCIVAollWoUTypDhmSsunv03SkEJCKVQtVf1UhTtYhIkSipVKMmksqEOROYMGdCOwckIpVC1V/VqImkMn/F/HYORkQqSckmFTN7BVgHbAW2uHttuhFVkG7doEsXVX+JSOJKNqlEPu3u+uZLWjOTSoqItJbaVKqVkoqIFEEpJxUHHjazRWY2LnulmY0zs4VmtnB1nnuDSBPyJJX+PfvTv2f/FAISkUpQytVfI9x9pZl9BPizmb3g7o81rHT3qcBUgNraWk8ryLJVUwOvvrrN4uknTU8hGBGpFCVbUnH3ldHfN4B7gAPTjajCqPpLRIqgJJOKmXU3sx4Nz4GjgGfSjarC1NTA2rWwZUvG4vEPjWf8Q+PTiUlEyl6pVn/tDNxjZhBi/K27P5RuSBWmYazKmjWw884fLq77T1068YhIRSjJpOLuLwGfSDuOihYfABlLKiIibVGS1V/SDjT/l7SnGTNg8GDo0CH8nTEj7YikSJRUqlXfvuGvkoq0VXMJY8YMGDcu9DZ0D3/HjVNiqVBKKtUqT0llSJ8hDOkzJMcOIjk0lzDc4ZJLYP36zP3Wr4dJk9o/Xik6JZVq1Se6GVdWUpl63FSmHjc1hYCkZDVVEpk0KXfC+OpXYf/9YccdYcWK3MddtqxYEUuKSrKhXtpB587Qo4eqv6RpDSWRhsTx6qtw1llwzz2w/fY5B9ACsGED7LorjBgRjrF27bbbDBxYvLglNSqpVLMcAyDHzR7HuNnbzIojlSxfSeS112D8+G1LIps2wd13w/z5YbbrXAYNggcegOuvD49u3TLXd+sGkycn/EakFCipVLMcSWXJmiUsWbMkpYCk3eVqE/nSl0I38912y1+SNYOXX4Zbbmk+YYwZA1OnNla59usXXo8ZU5z3JKlSUqlmmqql8jXVHvLWW7lLIlu3wrp1cO21sMsuuY/bUHXVkDAGDQqJZtCg3AljzBhYuDA8v/RSJZQKpqRSzZRUKluuUshXvgInnwwHHRS6lef79//gAzj/fLjmmsJKIq+8AvX14W++hDFoUCj9/O1vSbw7KVFKKtVMSaWy5eqZ9cEH8PvfhyTz3e/mn02hpSWRQpiFhnsllYqmpFLNamrgvffCF01k2C7DGLbLsPRikpbJVb21bFloHM/XM8sMFiyAyy4LVVxJlUQKMWIELF+u7sQVTEmlmsUnlYxMGTWFKaOmpBOPZGrNSPUzzwyliW99CzrlGTEQ78qbZEmkECNGhL8qrVQsJZVqpvm/0tVU0siVMM4+G849F666Cs47LwwwzK7ecofevWHxYpg2rbCuvEmWRJqz775hfJSSSsXS4MdqliOpnPH7MwDdAbLocg0q/OpX4dlnQ2P2hAnbJoyNG+EXvwjPe/cOAwxzeecdGDIkPCC0rSxbFkookyen2/OqY0c4+GAllQqmpFLNciSVFe/mmVJDkpWrEX3DBrjyyqb3Mwvdfbt3D6WbXO0m2dVbpdZ9d8QI+N//DaPsd9wx7WgkYQVVf5nZTjke2xU7OCkyVX+1v/p6eOihphvRV67MP4XJwIEhoUAodZTjSPURI0I13fz5aUdSvaKq1+EwPOlDF9qm8iSwGlgC/Dt6/oqZPWlmiQcl7WSnncJfJZXiiLeZDBgAZ5wBQ4fCZz8bluUycGCYM+uKKwofqd5ejexJOfDA0IlAVWDpiLfXFUGhSeXPwDHuXuPufYDPAvcDXwduLEpkUnydOoXqByWV5GU3tK9YEZaZhb+/+lXTSaMlI9Xbq5E9Kd26wfDhSirF1FQnkIkTt616TVChbSr/5e5fbXjh7g+b2TXufo6ZdS5GYGY2Cvgp0BG4xd2vKsZ5ql5NDaxe/eHLg/oflGIwFeLdd0OX3lz/cTdtgi98ITzv1KnpRvRSbA9JyogR8POfh84HnYvyFVK9cnUCOfts+N3vQieOIo8RKrSkssrMLjazQdHjIuB1M+sI1CcdVHTcGwglor2A081sr6TPI2wzqv7KI6/kyiObaSyW3L8E588P08L36xfm1col/h+6HEsZSRkxIiSUhvnApHBNlUI2b4YLL8zdc/APfwidPHr0KGp4hSaVLwD9gXujx8BoWUdgdBHiOhBY6u4vufsm4E7g+CKcRzRVS8vlG3R48MEwa1YoiTQ3EWO1O+SQ8FdVYNtq6fils86Co44K87n17AmrVuU+rllI4jfdtG3Va4IKSiru/qa7f9Pd948e57n7anff5O5LixDXbsDy2OsV0bIPmdk4M1toZgtXx6pvpIWyksrJs07m5FknpxhQGchVJ+0eOj6sWgW//GVhEzFWs759Q6cFJZVM+W7PfN11cP/9YdBrrvvbzJkTbpr29a839urMlms+tyIotEtxXzO72sweNLO/NDyKElGB3H2qu9e6e23fvn3TDKW8NSQVdwDWrF/DmvVrmtmpCuT6tbhuXRh8mK9Oeu3axqqFcu2Z1Z5GjIC//z1U/1WTpkoil1yS+/bM558Pxx0Hb7+d/7iPPhrmcpsypeD53BbBora9mW0VWv01A3gB2B34PvAK8ETSwcSsBAbEXvePlknSampCfev776cdSenI9Wtx7Njw6/rcc2G7PEO0squ2qrnNpBCHHhoS8fPPpx1JclozX9uXvwy1teGHx4omBh/Pmwf9++del+Z8blkKTSp93P1WYLO7P+ruZwGfKWJcTwB7mtnuZrY9cBpwXxHPV70aSnlqV2mUa7T7li2ht9aCBXDbbaraSkKlTS6Zr+rq17+Gp56C6dND9VT2Z2vz5rD+4IPD9Du5DBoU2kyuuqr05nPLUmhS2Rz9XWVm/21m+wM7FSkm3H0LcB7wJ+B5YJa7P1us81U1jarPtGxZ/kFh69fDpz6V+i/BirHHHqFDQzkllaZKIrl+jKxfH0q5w4aFzhzvvpv7uFu3wh13hG7WSYxfSlGh41QuN7NewHeA64GewPhiBQXg7g8CDxbzHMI2SeWI3Y9IMZh2NGNG4xiRAQPg1FPhhRfggQfy71Pqc2qVm3K7aVeu8R9f+Qo88khjySSfO+4IMzQfc0zuNrl4IzqU9/gld2/VAxjf2n2TfgwfPtyllZYscQf3229PO5L2M326e7du4X3HHz16uE+c6D5lyrbru3UL+0mypkwJ13f58rQjCaZPdx80yN0s/I3/m/fvv+1npuHRt69716651w0alHn8EvpsAQs94e/jtiSVZUkH09qHkkobvPVW+Bhcd13akSQn3xfDli3uCxa49+qV+z//wIHNH0OStXBhuPZ33pl2JLm/8Lff3v2II9z32Sd/QjFzr68vPGGU0Ger1JLK8qSDae1DSaUNtm5179jRfdIkd3cfNX2Uj5o+KuWg2iDfF0NtrXvv3vm/GBq+HKR9bd7s3r27+3nnpR1J/pKImfvRR+f/MZJdEimRhFGIYiSVttz50dta9SYloEMH6NPnwzaVDZs3sGFznps/lYNcjaWbNsE//wknnRTqtgvplinto1On0KupPdpVcjWyP/ss/OhHoXtzU915H3oIbrih4PEf1dyNvMmGejNbR+7kYUDXokQk7a+SpmrJNzCxvh5uvTU837o1s8EV1CU4TSNGwA9+ECY77NWrOOfI1ch+5pkfDvrlk58M537nnW33bUkjujRdUnH3Hu7eM8ejh7vrrpGVolKSyqxZoUdRLiU0OEyyjBgRkv6CBW07Tq6SyObNYdDgN76Rf2qdFStg0SKVRBLSluovqRTlnlTWrg3/uU89NXyZdOmSub7EBodJlk99Kty7vi1VYLkGHn7xi7DDDmHyylwlEAifnd2iaQX1YyMRSiqSkVSOHXIsxw45NuWAmhH/RbrzzmEQ3cyZ8P3vw+LFcMst+mIoJzvsAPvv33RSyTfo0B2eew6++c1tSyL19WGSxbvuKrwdTT822i7plv80Hur91UYTJ4YeYPX1aUfSvFy9u8zcf/CDtCOTthg/3r1LF/eNG7ddl69H30EHue+8c2E9+kpsfEipoMR6f0mlqKkJjdf5qghKSb5p5xsa4aU81dfDBx+Eqsvs6U8uvjh3j75//AOOPDL82++6a+7j5pruXSXYolJjuzRO1bJ6NSPvPQGAuWPnphZOXnV1+Xt3FfkWqVJEM2aEe9BA5o2nbrstNKKvzDNBuXuYpBHCLYmb69FX6tObVAiVVKT0J5VctSrcY/uTnwx16rlojEn5mjQJNmSNjdq0Cf7yF/jYx2DHHXPvpx59JUlJRUovqcQbZXv3Ds9vvz3cqCjXrVA1xqS8NVXKvP9+uP76kp/uXRqp+ksyk0qeYR7tJnuQ2jvvhO6mV18N3/52WNa9uwagVZKBA3PP8KtBh2VJJRUprZJKrob4rVvhpz9tfK1fpJVl8mQNOqwgKqlIGCew/fbw5puMHjk6vTjc1RBfjVQSqShKKhIaNqMBkF8/4EfpxOAeBrDlo4b4yqaeWRVD1V8SREll/eb1rN+8vvntk1RfH+ZmuuGGcGc8NcSLlK2SSypmdpmZrTSzuuhxTNoxVYUoqRwz4xiOmdGOl7y+Hs49N/Tquvji0NtHXUNFylapVn9d5+7XpB1EVampCYML6dd+56yvh3POCXN1TZwIl18eEomqQkTKVsmVVCQlffu2T++v+BiUnj1DQrn00saEIiJlrVSTynlm9rSZ/crMcg6nNbNxZrbQzBauXr26veOrPDU1sHYtHbcW8Yae2dOTv/8+bLcdfPzjSigiFSKVpGJmc8zsmRyP44GbgI8Cw4BVwLW5juHuU9291t1r+/bt237BV6qaGnBnh/c3F+8cuW71u3lzWC4iFSGVNhV3P7KQ7czsl8D9RQ5H4MMBkGP7H8s7e+SZ8bWtNAZFpOKVXEO9mfVz91XRyxOBZ9KMp2pESeXEvofCsEOTP358KvNsGoMiUjFKsU3lx2b2LzN7Gvg08O20A6oKUVJ5d+VLvLk+wQb7jRvDGJQzzoAhQ6Br18z1GoMiUlFKLqm4+5nuvq+77+fun4uVWqSYoqRy85+u4JRZp7T+OPHeXf37w157wY03wgUXwL/+Fe6boTEoIhWr5Kq/JCV9+gDQa10bGuqzZxhuuLnS//xPmGUYNAZFpMKVXElFUtK1K3TvTq/32pBUcvXuArj33tYfU0TKipKKNKqpaVtJRb27RKqekoo0qqmhd2uTysqVYSBjLurdJVI1lFSkUU0Ne/qOnFt7bsv2W7QIDjwwNL537py5Tr27RKqKkoo0qqmhz/v1nLrPqYXv87vfwaGHhlLKE0/Arbeqd5dIFVNSkUY1NdS/uZrl7yzPvT7eXXjQIDjlFBg9GvbfHx5/HPbdV7d9FalySirSqKaGDuve46xZORJB9mSQy5bB3XfDIYfAI4/ARz7S/vGKSMlRUpFG0QDInrm6FefrLrx8OXTpUuTARKRcKKlIoyip5Byrkq9b8PI8VWUiUpWUVKRRQ1LJ1a24X547Qqq7sIjEKKlIo3wllaefhnXrtt1e3YVFJIuSijSKkspJfQ9vXLZgARx+OPTqBT/+sboLi0iTNKGkNIomldyvY1TV9Ze/wOc+B7vsEnp4DRoEF16YYoAiUupUUpFGs2bhZvj3vhe6CB91FOy+O/z1ryGhiIg0Q0lFgmgcirljAKtXh/Eo3/xm/kZ6EZEsSioS5BqHUl8PV1yRTjwiUpaUVCTQtPUikoBUkoqZfd7MnjWzejOrzVo3wcyWmtliMzs6jfiqUr7xJhqHIiItkFZJ5RngJOCx+EIz2ws4DdgbGAXcaGYd2z+8KjR5chh3EqdxKCLSQqkkFXd/3t0X51h1PHCnu29095eBpcCB7RtdlRozBqZOZcOuH8E1DkVEWqnUxqnsBiyIvV4RLduGmY0DxgEMVBVNMsaMoauSiIi0QdGSipnNAXbJsWqSu/+hrcd396nAVIDa2lpv6/EkqPtPHQDDdhmWahwiUp6KllTc/chW7LYSGBB73T9aJu1k/EPjAZg7dm6qcYhIeSq1LsX3AaeZWWcz2x3YE3g85ZhERKRAaXUpPtHMVgAHAQ+Y2Z8A3P1ZYBbwHPAQ8A1335pGjCIi0nKpNNS7+z3APXnWTQbUj1VEpAyVWvWXiIiUsVLrUiwpu+IIzfUlIq2npCIZDh5wcNohiEgZU/WXZJi3fB7zls9LOwwRKVMqqUiGiY9MBDRORURaRyUVERFJjJKKiIgkRklFREQSo6QiIiKJUUO9ZJgyakraIYhIGVNSkQya8l5E2kLVX5JhzktzmPPSnLTDEJEypZKKZLj8scsBOHKP1twOR0SqnUoqIiKSGCUVERFJjJKKiIgkRklFREQSo4Z6yXDzsTenHYKIlLG07lH/eTN71szqzaw2tnywmW0ws7ro8Ys04qtmQ2uGMrRmaNphiEiZSquk8gxwEpDrZ/GL7j6sfcORBrMXzwbguKHHpRyJiJSjVJKKuz8PYGZpnF6acO38awElFRFpnVJsqN/dzP5pZo+a2aFpByMiIoUrWknFzOYAu+RYNcnd/5Bnt1XAQHdfY2bDgXvNbG93fzfH8ccB4wAGDhyYVNgiItIGRUsq7t7ieT7cfSOwMXq+yMxeBIYAC3NsOxWYClBbW+tti1ZERJJQUtVfZtbXzDpGz/cA9gReSjcqEREpVCoN9WZ2InA90Bd4wMzq3P1o4DDgB2a2GagHvubub6URY7W6/cTb0w5BRMpYWr2/7gHuybH8buDu9o9IGgzoNSDtEESkjJVU9Zekb+YzM5n5zMy0wxCRMqVpWiTDTQtvAuDUfU5NORIRKUcqqYiISGKUVEREJDFKKiIikhglFRERSYwa6iXDXaPvSjsEESljSiqSoaZbTdohiEgZU/WXZJhWN41pddPSDkNEypSSimRQUhGRtlBSERGRxCipiIhIYpRUREQkMUoqIiKSGHUplgwPjnkw7RBEpIwpqUiGbtt1SzsEESljqv6SDDc+cSM3PnFj2mGISJlSUpEMs56dxaxnZ6UdhoiUKSUVERFJTCpJxcyuNrMXzOxpM7vHzHrH1k0ws6VmttjMjk4jPhERaZ20Sip/BvZx9/2AJcAEADPbCzgN2BsYBdxoZh1TilFERFoolaTi7g+7+5bo5QKgf/T8eOBOd9/o7i8DS4ED04hRRERarhS6FJ8FzIye70ZIMg1WRMu2YWbjgHHRy41m9kzRIkxODfBm2kEUoMa+bGURJ2VyPVGcSSmHGKF84hya9AGLllTMbA6wS45Vk9z9D9E2k4AtwIyWHt/dpwJTo+MsdPfaNoTbLhRnshRnssohznKIEcorzqSPWbSk4u5HNrXezMYCxwJHuLtHi1cCA2Kb9Y+WiYhIGUir99co4CLgc+6+PrbqPuA0M+tsZrsDewKPpxGjiIi0XFptKj8HOgN/NjOABe7+NXd/1sxmAc8RqsW+4e5bCzje1OKFmijFmSzFmaxyiLMcYoQqjtMaa55ERETaRiPqRUQkMUoqIiKSmJJMKmY2KpqmZamZXZJjfWczmxmt/4eZDY6tyznNS3PHbM84zez/mtkiM/tX9PczsX3mRsesix4fSSnGwWa2IRbHL2L7DI9iX2pmP7OoYSylOMfEYqwzs3ozGxatS/RaFhjnYWb2pJltMbNTstZ9ycz+HT2+FFuexvXMGaeZDTOz+Wb2rIVplE6NrZtmZi/HruewtOKM1m2NxXJfbPnu0WdkafSZ2T6tOM3s01mfzw/M7IRoXRrX83wzey76t33EzAbF1iXz+XT3knoAHYEXgT2A7YGngL2ytvk68Ivo+WnAzOj5XtH2nYHdo+N0LOSY7Rzn/sCu0fN9gJWxfeYCtSVwLQcDz+Q57uPAfwEG/BH4bFpxZm2zL/BiMa5lC+IcDOwH/AY4JbZ8J+Cl6O+O0fMdU7ye+eIcAuwZPd8VWAX0jl5Pi2+b5vWM1r2X57izgNOi578Azk0zzqzPwFtAtxSv56dj5z+Xxv/viX0+S7GkciCw1N1fcvdNwJ2E6Vvijgd+HT2/Czgiyp75pnkp5JjtFqe7/9PdX4uWPwt0NbPObYwn0RjzHdDM+gE93X2Bh0/cb4ATSiTO06N9i6XZON39FXd/GqjP2vdo4M/u/pa7ryXMfzcqreuZL053X+Lu/46evwa8AfRtYzyJx5lP9Jn4DOEzAuEzc0KJxHkK8EfPHEaRpELi/H+x88enyErs81mKSWU3YHnsda6pWj7cxsMcYu8AfZrYt5BjtmeccScDT7r7xtiy26Li8HfbWBXS1hh3N7N/mtmjZnZobPsVzRyzveNscCpwR9aypK5loXG2dN+0rmezzOxAwi/eF2OLJ0dVJ9cl8EOorXF2MbOFZragoUqJ8Jl42xvnFiyZ60koYWd/PtO8nmcTSh5N7dviz2cpJpWqYWZ7Az8CzoktHuPu+wKHRo8z04iNUO0x0N33B84HfmtmPVOKpVlm9ilgvbvH54ArlWtZdqJfqLcDX3b3hl/fE4CPAwcQqkkuTim8BoM8TIXyBWCKmX005Xjyiq7nvsCfYotTu55mdgZQC1yd9LFLMakUMlXLh9uYWSegF7CmiX2LMf1LW+LEzPoD9wBfdPcPfwm6+8ro7zrgt7RtluZWxxhVIa6JYllE+LU6JNq+f2z/1K9lZJtfgQlfy0LjbOm+aV3PvKIfDw8Q5un7cIJXd1/lwUbgNtK9nvF/35cI7Wf7Ez4TvaPPSIuPWYw4I6OBe9x9c8OCtK6nmR0JTCLMaLKxmX1b/vlMqpEoqQdhlP9LhIb2hsamvbO2+QaZjbazoud7k9lQ/xKh8arZY7ZznL2j7U/Kccya6Pl2hHrhr6UUY1+gY/R8j+iDtJPnbrg7Jq1rGb3uEMW3R7GuZaFxxradxrYN9S8TGkF3jJ6ndj2biHN74BFgfI5t+0V/DZgCXJVinDsCnaPnNcC/iRqlgd+R2VD/9bTijC1fAHw67etJSLwvEnXGKMbns9VvoJgP4BjCzbteJPxaAvgBIbMCdIk+OEujNxz/MpkU7beYWC+FXMdMK07gUuB9oC72+AjQHVgEPE1owP8p0Rd7CjGeHMVQBzwJHBc7Zi3wTHTMnxPNzJDiv/lIwlQ/8eMlfi0LjPMAQr3z+4Rfzc/G9j0rin8poVopzeuZM07gDGBz1mdzWLTuL8C/olinAzukGOfBUSxPRX/Pjh1zj+gzsjT6zHRO+d99MOFHT4esY6ZxPecAr8f+be9L+vOpaVpERCQxpdimIiIiZUpJRUREEqOkIiIiiVFSERGRxCipiIhIYpRUpKJkzVxbZ7EZrMudme1vZrdGz8ea2c+z1s81s9om9r/TzPYsdpxS3dK6nbBIsWxw92G5VkRzf5k3TjtSbiYCl7dh/5uAi4CvJhOOyLZUUpGKZuG+MIvN7DeEAVwDzOxCM3simsjv+7FtJ5nZEjP7m5ndYWYXRMs/LAGYWY2ZvRI972hmV8eOdU60fGS0z11m9oKZzWiYzNLMDjCzeWb2lJk9bmY9zOyx+L00ovN/Iut99AD2c/enCnjPn4uV1Bab2cvRqr8CR8amMBFJnD5cUmm6mlld9Pxl4NvAnsCX3H2BmR0VvT6QMO3EfWZ2GGEk9GnAMML/iycJI/KbcjbwjrsfEM0w+3czezhatz9h2qDXgL8Dh5jZ48BM4FR3fyKaY2sDcCswFhhvZkOALjmSR8Oo5rhTzWxE7PXHANz9PuA+ADObBTwaLa83s6XAJwp4byKtoqQilSaj+itqU3nVGydGPCp6/DN6vQMhyfQgTPi3PtrvwzsJNuEoYD9rvNNfr+hYm4DH3X1FdKw6wlQd7wCr3P0JAHd/N1r/O+C7ZnYhYaqMaTnO1Q9YnbVsprufF3uvc+MrzewiwvW4Ibb4DcLNt5RUpCiUVKQavB97bsCV7n5zfAMzG9/E/ltorCruknWsb7p7fDpzzGwkEL8/zlaa+L/m7uvN7M+EGyqNBobn2GxD1rmbFM1E+3ngsKxVXaJjiRSF2lSk2vwJOMvMdgAws90s3Lv+MeAEM+satV8cF9vnFRq/6E/JOta5ZrZddKwhZta9iXMvBvqZ2QHR9j1i7Ru3AD8DnvBw571szxNVbzXHwn3HbwA+7+7ZCWQI21ajiSRGJRWpKu7+sJn9H2B+1Hb+HnCGuz9pZjMJs96+ATwR2+0aYJaZjSPcZ6TBLYRqrSejhvjVNHGrVXffZGanAtebWVdCieFIwr3WF5nZu4T7auTa9wUz62VmPTzcH6YpYwl3QLw3eo+vufsxZrYzoTrsP83sL9JqmqVYJAczu4zwZX9NO51vV8KNpj6er8uzmX0bWOfut7TyHN8G3nX3W1sdqEgzVP0lkjIz+yLwD8L9L5oaQ3MTmW01LfU28Os27C/SLJVUREQkMSqpiIhIYpRUREQkMUoqIiKSGCUVERFJjJKKiIgk5v8D8ZVR1RH7rk8AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "cross = cross.rebin(0.0050)\n", "lag, lag_e = cross.time_lag()\n", "plt.figure()\n", "\n", "# Plot lag-frequency spectrum.\n", "plt.errorbar(cross.freq, lag, yerr=lag_e, color='r', fmt=\"o-\")\n", "\n", "# Find cutoff points\n", "v_cutoff = 1.0/(2*delay)\n", "h_cutoff = lag[int((v_cutoff-cross.df*2)*1/cross.df)]\n", "\n", "plt.axvline(v_cutoff, color='g',linestyle='--')\n", "plt.axhline(h_cutoff, color='g', linestyle='-.')\n", "\n", "# Define axis\n", "plt.axis([0,0.2,-20,20])\n", "plt.xlabel('Frequency (Hz)')\n", "plt.ylabel('Lag')\n", "plt.title('Lag-frequency Spectrum')\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Energy Dependent Impulse Responses" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In practical situations, different channels may have different impulse responses and hence, would react differently to incoming light curves. To account for this, stingray an option to simulate light curves and add them to corresponding energy channels.\n", "\n", "Below, we analyse the lag-frequency spectrum in such cases. \n", "\n", "We define two delta impulse responses with same intensity but varying positions, each applicable on different energy channels (say '3.5-4.5 keV' and '4.5-5.5 keV' energy ranges). " ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "delays = [10,20]\n", "h1 = sim.simple_ir(start=delays[0], width=1)\n", "h2 = sim.simple_ir(start=delays[1], width=1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, we create two energy channels to simulate light curves for these two impulse responses." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [], "source": [ "sim.simulate_channel('3.5-4.5', var, h1)\n", "sim.simulate_channel('4.5-5.5', var, h2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compute cross-spectrum for each channel." ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [], "source": [ "cross = [Crossspectrum(lc, var).rebin(0.005) for lc in sim.get_channels(['3.5-4.5', '4.5-5.5'])]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Calculate lags." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "lags = [c.time_lag() for c in cross]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Get cut-off points." ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [], "source": [ "v_cuts = [1.0/(2*d) for d in delays]\n", "h_cuts = [lag[int((v_cutoff-0.005*6)*1/0.005)] for lag, v_cut in zip(lags, v_cuts)]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Plot lag-frequency spectrums." ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZUAAAEWCAYAAACufwpNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABZAUlEQVR4nO2dd3gU1frHP286kNBCCdJRQkeaVEVaAAVEpHv1gohYfurFfu0oYrkiKioKNlCQIohSlSYiTYoEpYcSJHRCCyUJIef3x5mFTUjZJLs7m+R8nmeenZ1z5px3zs7Od057jyilMBgMBoPBHfjZbYDBYDAYCg5GVAwGg8HgNoyoGAwGg8FtGFExGAwGg9swomIwGAwGt2FExWAwGAxuw4iKId8iIiNEZLLddhQmTJkbssOIiocQkVgRuSgi55y2j+22KztERInIecveeBFZKiL97bbLk4hINeu6A7KIM0JELqX7PZ/1pp0FHRG5X0R2iEiCiBwVkQUiEubB/LL93Q05xxSmZ+mhlFriyQxEJEApleLmZG9USu0WkTLAbcDHIlJbKfWam/PJb0xXSt2TVQQR8VdKXfaWQQUFEbkVeBPoqpTaJCKlgR42m+Wp/1eBxtRUbEBEBovIShEZLSKnRGSfiNzmFF5CRL4UkcMiclBE3hARf6dzV4nI+yISD4wQkXARmSsiZ0VkvRV/pRX/ExF5L13+c0TkiezsVEqdUEp9CzwMPC8i4Tmw72MROWO9eXbMwbVlVS7VReQ36012MVAm3XW1FJHVInJaRDaLSDunsOUiMtKyLUFEFlmiCbDC+jxt1UBaZVc2TulOFJFPrbfq80B7EblORGaJyHHrGh53il/EOueUiGwTkWdEJM4pXInIDenSf8Ppe3cRibaucbWINHQKixWRp0XkL6vsp4tIiFN4T+vcsyKyR0S6ikhfEdmY7pqeFJGfXLz+70XkiJXfChGp5xSW6X2ZATcBa5RSmwCUUieVUpOUUglO5fCZiCy2fr/fRKSqU161rbCTIrJTRPqlK/P3RGS/ZedKESlCBr+7ZPz/StPkJ+lqONa99Yb1e5yzrjlcRKY4XXs1V8qzQKCUMpsHNiAW6JRJ2GDgEvAA4I9+aB8CxAqfDYwHigHlgHXAg07npgCPoWuaRYBp1lYUqAscAFZa8ZtbaftZ38sAF4DymdimgBvSHQu08rwtB/Y9YZ3XHzgDlHbx3KzKZQ0wBggG2gIJwGQrrCIQD9yOflmKsr6XtcKXA3uASKvMlgNvW2HVrOsOyOL3HOHIK93xidb1tbHyLQpsBF4BgoAawF6gixX/beB3oDRQGdgCxGVW/lb6b1j7jYFjQAurfAah77Ngp3tuHXCdlf524CGn++CMVS5+VnnVtsryJFDHKc9NQG9XygEYAoRZ6XwARDuFZXpfZpDuLcBF4DWrLIMzKOcE63cPBj7k6j1ezEr7PvR/ojFwAqhrhX9i/d4VrXJrbaVxze9Oxv+v9Nec5jwr7d3A9UAJYBuwC+hkpfEN8LXdzySvPfvsNqCgbtYf/Bxw2ml7wAobDOx2ilvUukkjgPJAElDEKXwg8KvTuf84hfmjH8S1nI694fznRT9coqz9R4EFWdh9jahYx48A/3LRvitCYB1bB9zr4rmZlUsV689ezCn8O66KynPAt+ls/gUYZO0vB15yCnsE+Nnav+bhksH1jwCS0/2e16Efdt84xWvh/PtYx553PFTQAtPVKWwYrovKp8DIdGnvBG51uufucQr7H/CZtT8eeD+Ta/sUGGXt1wNOke6hnq4crhFXK6ykZX8JXLgvMzj/NmCuVbbn0C8Q/k7lMM0pbihwGS3M/YHf06U1HngVLaAX0U266fO75ncn3f8ro2tOf551b73oFP4esNDpew+cxLagb6ZPxbPcqTLvUzni2FFKXRAR0H+U0ug3/MPWMdB/jANO5zrvl0W/DWUWDjAJuAdYbH1+mJOLEJFAK5+TQFUX7DuorH+TxX70A9iVczMrlzLAKaXU+XTpVrb2qwJ9RcS5HT4Q+DWjtNG1tdBMLzpjZqh0fSqWfc72VwWuE5HTTsf80bUT0OXgHH9/DvKvCgwSkcecjgVZaTpIf42OsMrAgkzSnQRMFZGX0OI/QymVJCL/Qj+cQT+0b3M+SXSz5SigL/r+SLWCyqDf8LO7L9OglFoILBQRP6A98D1aNB02HHCKe05ETnL1vmqRrswDgG8tW0LQtVRXydLOTDjqtH8xg+85vdfyLUZUfI8D6Lf5MirzDkLnB/Zx9Bt8JXSVG64+aB1MBraIyI1AHeDHHNrU08pjHfohlp19FUVEnISlCjAH164tMw4DpUSkmJOwVOFqWRxA11QeyGG6kLY8c4Pz+QeAfUqpmpnEPYz+fbZa36ukC7+ArqE5iAAcfS4H0DWKUbmw8QC6eeYalFJrRSQZ3QR1t7WhlJoCTMkizbvR90YndC2pBLqWI7h2X2aIUioVWCoiy4D6TkFXzhcRxwvYIevaflNKRaVPyxKoRPS1b06fVWYmpPt+nmt/E0MmmI56H0MpdRhYBLwnIsVFxE9Erhc9Oiaj+JeBH9AdikVFpDbw73Rx4oD16De3WUqpi67YIiKlrbfVT4B3lFLxLtpXDnhcRAJFpC9ayBbk9NrSXcN+YAPwmogEicjNpB0dNBnoISJdRMRfREJEpJ2IVHLhUo+j37JruBA3O9YBCSLynNVB7C8i9UXkJit8BnrQQynLtsfSnR8N3G2d1xVwLpvPgYdEpIVoiolIN3Ft2O2XwH0i0tEq94rWveLgG+Bj4JJSKrPO9PSEoV8S4tEP3TcdAa7cl86IHkQwwCoXEZHm1rWvdYp2u4jcLCJBwEhgrVLqADAPiBSRe617LlBEbhKROpZAfQWMET2Awt/qkA/G9d89GmgrIlVEpAS6OdOQCUZUPMtcSTuvYbaL5/0bXSPYhn7zmwlUyCL+o+i3xCNo4ZiK/rM7MwloYIVnx2YROYfufBwKPKGUeiUH9v0B1ER3lo4C+iil4nN5bc7cje6zOIluL//GEWA9XHoCL6AfFgeAZ3DhHldKXbDsXCV6VFVLF+3JKK3LQHegEbAPXQZfoH8f0B3R+62wRVz7e/wHLZan0X1YPzqlvQE9iOFjdNntRvcBuGLXOnRH9vvoDvvf0M1GDr5F1wpyMrHxG+taDqJ/z7Xpwl25Lx2cQl9bDHDWsuNdq7bk4Dv0734SaIpuykXpEWKdgQHomssR4B10ZzzA08Df6Berk1aYn6u/u1JqMTAd+As9CGNeFmVS6HGMqjEUIETkHSBCKTXI6Vhb9B+1qvLgjy4ig4GhSqmbPZVHQUL0sOfJSilXalSetKMIemRZE6VUjIfyuOa+zMG5E9EDGl5yu2EGt2JqKgUA0WP0Gzo1G9yPHrrrCA9EvwF/4UlBMeRrHgbWu1NQsrsvDQUT2zrqRaQyuvpcHt0xNkEp9aGIjEBXg49bUV9QSmU2asWgCUM3LVyHHnXyHvATgIjUQfdFbEY3fxgMaRCRWHTn+p1uTjrT+9JQcLGt+UtEKgAVlFJ/Wh2NG9E3dT/gnFJqtC2GGQwGgyHX2FZTsUYCHbb2E0RkO3rGq8FgMBjyKT7RUS/aL84K9OiTJ9EjWs6im22eUkqdyuCcYejZyBQrVqxp7dq100e5wtmks8TEx1CrTC1Cg+yZg7QzficAtcJr2ZK/y+zUdlLLx+00GAx5ZuPGjSeUUmXdmabtomJNYvoNPanrBxEpjx6GqdBj0SsopYZklUazZs3Uhg0bMg3fc3IPN3x0A1/e8SVDGmeZlMdoN7EdAMsHL7clf5dp105/Ll9upxUGg8ELiMhGpVQzd6Zp6+gva1TSLGCKUuoHAKXUUaXUZWvS0udoR3h5omrJqgT6BRIT75GRkgaDwWCwsHP0l6Bn+W5XSo1xOl7B6m8B6IX24ponAvwCqFGqBjEn7ROVl9rmk+H1L+UTOw0Gg09ip++vNmjndX+LSLR17AVgoIg0Qjd/xQIPuiOzmuE12RW/K/uIHqJTjU625Z0jOuUTOw0Gg09i5+ivleix8enxyJyUyNKRLN27lFSVip94v9Uv+kg0AI0iGnk97xwRHa0/GzWy0wqDwZBPKTReimuG1+RiykUOJRyiUnHve8QY/vNwIB901A8frj9NR73BYMgFhcZNS83S2hO5nU1gBoPBUNApNKISGR4JYEaAGQwGgwcpNKJSsXhFQgJCbB0BZjAYDAWdQiMqfuLHDaVvMM1fBoPB4EEKTUc96Cawbce32ZL3mx3fzD6SL/BmPrHTYDD4JIVKVGqWrsncnXO5nHoZfz9/r+bdunJrr+aXa1rnEzsNBoNPUmiav0CLyqXUS+w/s9/rea8+sJrVB1Z7Pd8cs3q13gwGgyEXFKqaivMIsBqlang17xeWvgDkg3kqL2g7zTwVg8GQGwpXTSVcz1UxI8AMBoPBMxQqUSlfrDyhQaFmBJjBYDB4iELV/CUiRIZH5qmmcuHSBdYdXMfqA6spFliMqOujqFOmDtrpssFgMBRuCpWogO6s33Ao8wW90hN/IZ5VB1bx+/7fWXlgJRsPbeRS6qU0ca4Lu46oGlFE1YiiU41OlA8t726zDQaDIV9g+8qP7iCsephq+mrTLON0j+zO062f5uVlL/PG72/weY/PGdpkKCcunKDPjD5p4iZfTmb/6f2cTjzNhZQLAPiLPy0rtaRphaYsi13GCze/QOvKrfn2r28Z+8dYTiWeIiU1BYBigcUoVaQUZYuWpXhwcQB61e7FrdVupUhAER6c9yBvdnyT1pVbs/rA6iud+FmRPv747uOpVaYWc3fO5b0172V7fvr4M/vNpEzRMkyMnsjE6IlXI547pz9D0y67nD6+Y8DB6NWjmbdrXrb5O8dfE7eGWf1mAfD8kudZE7cmy3PDi4aniR9/MZ4JPSYAMGzusGybMyPDI9PEDy8Szlud3gKg94zexF+Iz/L8VpVapYnfqlIrnm79NHB1Rc+scNx7jviDGw1mcKPBGd57GZE+/lOtnqJHrR7sPLGTB+dlvzJE+vg+e+9lgrn3PHfv/Xbfb25f+bHQ1VQcI8COnz+eYfix88eIiY/hsrpMqZBSlA8tT4ngEjxy0yM80PQBTlw4weajmwkNCqVqyar0rduXJXuXAJCQnMCpi6c4lXiKg2cPEnc2jluq3IKf+FGjVA0aRTRi54mdXrvWXJFOTAwGgyEn2FZTEZHKwDdAefSCXBOUUh+KSGlgOlANvUhXP6XUqazSym6NemfWxq2l1ZetmDtwLt0ju185Hn8hnkcWPMKMrTNoWaklk+6cdEWAcsO49eP4vwX/x+GnDhMRGnFFeHx+sa4l2k6zWJfBUPDxxBr1dtZUUoCnlFJ/ikgYsFFEFgODgaVKqbdF5L/Af4Hn3JVpRi7w5++az9C5Q4m/EM+oDqN4ts2zBPjlrWjCgsIAOJesm5PeWPEGkA9E5Q1tpxEVg8GQG+xc+fEwcNjaTxCR7UBFoCfQzoo2CViOG0UlvGg4pYuUJiY+hrNJZ3nylyf5ctOXNCjXgJ//9TM3RtzolnxCg3QzkkNUDAaDoTDgE30qIlINaAz8AZS3BAfgCLp5zK3ULF2TZbHLuPGzG/nnzD/8t81/GdFuBMEBwW7LIyw4bU3FYDAYCgO2T34UkVBgFjBcKXXWOUzpDp8MO31EZJiIbBCRDcePZ9zpnhmR4ZHsit9FgF8Av9/3O291esutggKmpmIwGAonttZURCQQLShTlFI/WIePikgFpdRhEakAHMvoXKXUBGAC6I76nOQ7vOVw6pSpw+MtHqdYULE8XEHmGFExGAyFEdtERfQU9C+B7UqpMU5Bc4BBwNvW50/uzrtJhSY0qdDE3cmmIb2ojO8+3qP5uY3x+cROg8Hgk9hZU2kD3Av8LSLR1rEX0GIyQ0TuB/YD/ewxL284RCUhKQGAWmVq2WmO69TKJ3YaDAafxM7RXyuBzBxmdfSmLZ4gfU1l7s65APSo1cM2m1xirraTHj5up8Fg8El8YvRXQSTYPxh/8b8iKg53Fj4vKu9ZbjeMqBgMhlxg++ivgoqIEBoUajrqDQZDocKIigcJCw4zomIwGAoVRlQ8SGhQKOcuGVExGAyFByMqHsQ0fxkMhsKG6aj3IM6i8m2vb222xkW+zSd2GgwGn8SIigcJDQrl4NmDAFQuUdlma1ykcj6x02Aw+CSm+cuDONdUpm+ZzvQt0222yAWmT9ebwWAw5AJTU/EgoYFXReXTDZ8C0L9+fztNyp5PtZ3093E7DQaDT2JqKh7EdNQbDIbCRsEUFaUgMdFuK66Iil1LNhsMBoO3KXiikpwMVarAqFF2W0JYcBgKxcWUi3abYjAYDF6h4IlKUJAewbR4sd2WmDVVDAZDoaNgdtR36qRrKqdOQalStpnhLCoz+820zY4cMTOf2GkwGHySgldTAYiKgtRUWL7cVjOcRaVM0TKUKVrGVntcokwZvRkMBkMuKJii0rIlhIba3gTmvFDXxOiJTIyeaKs9LjFxot4MBoMhF9gqKiLylYgcE5EtTsdGiMhBEYm2tttznHBgILRr5zOici75nBEVg8FQKLC7pjIR6JrB8feVUo2sbUGuUu7UCXbvhtjYPJiXN0xHvcFgKGzYKipKqRXASY8kHhWlP22srRhRMRgMhQ27ayqZ8aiI/GU1j2U4fEtEhonIBhHZcPz48Wsj1KkD110HS5Z42tZMCQsKA4yoGAyGwoMvisqnwPVAI+Aw8F5GkZRSE5RSzZRSzcqWLXttBBFdW1m6VI8EswFTUzEYDIUNn5unopQ66tgXkc+BeblOrFMnmDQJNm2Cpk3dYV6OCAkIwU/8OJd8jgX/yl3XkNdZkE/sNBgMPonP1VREpILT117AlsziZkunTvrTpn4VEbni/6toYFGKBha1xY4cUbSo3gwGgyEX2D2keCqwBqglInEicj/wPxH5W0T+AtoDT+Q6g4gIaNDA1n6V0KBQEpITGLd+HOPWj7PNDpcZN05vBoPBkAtsbf5SSg3M4PCXbs0kKgo++QQuXoQiRdyatCs4aiozts4A4JGbHvG6DTlihraTR3zcToPB4JP4XPOX2+nUCZKS4PffbcnerKliMBgKEwVfVNq21Z6LbepXMaJiMBgKEwVfVIoVg9atbROVsKAwIyoGg6HQUPBFBXS/yubNcOyY17M2NRWDwVCYKByi4hhavHSp17N2iMrywctZPni51/PPMcuX275kgMFgyL8UDlFp2lQv1mVDE5ipqRgMhsJE4RAVf3/o0EGLilJezdoxT+XdVe8yevVor+adK0aP1pvBYDDkgsIhKqD7VeLiYNcur2YbGhRKqkpl7q65zNuVe48zXmPePL0ZDAZDLig8omKTyxaHU8nLqZe9mq/BYDDYQeERleuvh+rV7RMVZUTFYDAUfAqPqIBuAvv1V7h0yWtZOtZUMTUVg8FQGCh8opKQAOvWeS1LR00lwC+AIoHe9z2WY4oUscVHmsFgKBj43HoqHqV9e71415Il0KZNzs69eFF39B84oLe4OP3wrV1brzJZtSr4XavRDlF5vf3rRF0flbM8ExPh8mXtFcBbLFzovbwMBkOBo3CJSni4nrOyeDG8+mrm8f76CyZOhD17rorIiRNZpx0SArVqaZFxCE3duoSWDgFysPpjfDzMnw8//gi//KKFpW5daNHi6lavnh4mnRVK6VrZ4cNQoQIUL+5a/gaDwZAHRHl53oYnaBYWpjZkt7Jj9+7w9NPwwgvw1lt6zZCHH9Zi0aePjnPypK6BnDqlazRFi0JwsN5atYKePfXD+a234NlnoV07WLQIRoyACxeubomJV7LdEy7c8JjirpMRNKrRipevHwLvvqvTaN0aVq+GJ5/Udpw4AWfO6BODgqBMGQgMhLNnISVFiwToGlHTprrmlZwMc+dqT8zJyWk/nZdRrlkTbr1Vi9/atXrYcPnyWjwnTrwab/9+/Vm1atrymzlT2+OI75h1P3q0a0OQneOvWQOzZunvzz+vv2dFeHja+PHxMGGC/j5sWPbDxCMj08YPD9flD9C7t04vK1q1Shu/VSt9L4G+B7LDce854g8erDfney8r0sd/6ino0QN27oQHH8z+/PTx33zz6r33wgvZn58+/vjx+gVq7lx4L8PVvtOSPn76eyk7zL13Nb6b7z357beNSqlm2SfkOrbWVETkK6A7cEwpVd86VhqYDlQDYoF+SqlTbss0Kkr/SDt26O9JSXDkiBaT8+f1w7x6dbjuOghwKp4ePeDf/9Y/RpEiOl54ODRpomsCzqSm6uayc+cIvXwGOMzfHCE+ejYvD5+tReGhh6BlSz1wYPdufV7RolCliv4DhYWlTXPUKChXDiZPhi+/1ML1/vtXBx2IaPELCoLQUG1bUJDeLl7UC5bNmqUFE/RouJYt9bUkJFzNzxGeXlQMBoPBBWytqYhIW+Ac8I2TqPwPOKmUeltE/guUUko9l1U6zZo1Uxs2bHAt06Qk7bKlb1/dTDV2rBaVBg30G93AgfpB7CbOJ58n9K1QapSsTuXAMiwPeww2boQNG2DLFp3vnXfqWtANN+Qs8aQkiI2FsmX1NYlkHV8piInRNZU//tCfmzdrwYmL03Ecbz/G/5fBUOARkYJVU1FKrRCRaukO9wTaWfuTgOVAlqKSI4KD9Ror33yjv3fuDJMm6RpMdg/lXFA0sCiCcFml6prIvffqzR0EB+tmBVcR0dXxyEhd6wJ4+WXdvHH5cvb9NAaDwZANvthRX14pddjaPwKUzyiSiAwDhgFUqVIlZzk8/7x+sA4dCg0b5sHU7BERQoNCfXeeSkSEbq47cUL3sRgMBkMe8EVRuYJSSolIhu1zSqkJwATQzV85SvjWW/XmJUKDQgnwDyC8aLjX8nSZiAj9eeSIFpVwH7TRYDDkG3xRVI6KSAWl1GERqQB4f2UtNxMaFEqz65rxXe/v7DblWhy1kyNH4MYbr450MRgMhlzgizPq5wCDrP1BwE822uIWfHpNFeeaisFgMOQRW0VFRKYCa4BaIhInIvcDbwNRIhIDdLK+52tCg0LZfHQzzy953m5TriW9qDz/vN4MBoMhF9g9+mtgJkEdvWqIhwkNCuXkxZOsictmopUdhIZqNzBHj+rv2U0GMxgMhizwxeavAodPj/4CXVsxzV8Gg8ENGFHxAqFBob69nkr58kZUDAaDWzCi4gVMTcVgMBQWjKh4gbCgMFJVKpWKV7LblIyJiLjap1Kpkt4MBoMhF/jiPJUCR2hQKArFl3d8abcpGRMRoT00JyVph5UGg8GQS0xNxQs4Fury+bkqx/L9PFODwWAzRlS8gENUnln8jM2WZILzrPrhw/VmMBgMucA0f3kBh6j8fexvmy3JBOcJkNHRtppiMBjyN6am4gUcouKzI8AcouLorDcYDIZcYkTFC1wRFV+dq+Lc/OWDKKVYtm8ZW49ttdsUg8GQDUZUvIDP11SCg/XKkT4mKkop5u6cS/MvmtPxm44MmzfMbpMMBkM2mD4VLxAWrNd/L1usrM2WZIFjVn1kpN2WkKpS+WnHT4xcMZJNRzZRvWR1ml3XjB0ndthtmsFgyAZTU/ECjppK/3r9bbYkCxyz6idM0JsNpKpUZm6bSePxjblrxl0kJCfwdc+v2fnoTu6ufzcnL54k/kK8LbYZDAbXMDUVL+Dz81RAi8qGDbZkrZTi+23f8/pvr7P1+FYiwyP55s5vGNhgIAF++hatGV4TgJiTMb65gqbBYABMTcUrFA0sCsDULVNttiQLHDWVYcP05iWOnjtKz2k96T+zP6kqle/u+o5tj2zj3hvvvSIoAJHhulluV/wur9lmMBhyjs/WVEQkFkgALgMpSqlm9lqUe/zEDz/x49h5H56xXr48nDsH27eDv79Xsvxxx488MPcBEpISeL/L+zzW/DH8/TLOu3rJ6viLvxEVg8HH8VlRsWivlDphtxHuwF/8fXf0F1ydq5KcDEWKeDSrs0ln+c/P/2Fi9EQaRzRm8l2TqVu2bpbnBPoHUqNUDSMqBoOP4+uiUmDw9/P33XkqcFVULl3yqKis2L+Cf8/+NwfOHuDFW17klVtfIcg/yKVzI8MjjagYDD6OL/epKGCRiGwUkWsa+UVkmIhsEJENx48ft8G8nJGvaioeICkliWcWPUO7ie0I8Atg5X0reaPDGy4LCmhRiTkZg1LKIzYaDIa848s1lZuVUgdFpBywWER2KKVWOAKVUhOACQDNmjXz+adMiZAS+IkPa7hDVMqVg4YN3Zr09uPb6T+zP38f+5sHmz7I6M6jr4yIywk1S9fkwqULHEo4RMXiFd1qo8FgcA8+KypKqYPW5zERmQ00B1ZkfZbv0rB8Q05ePGm3GZlTpgyIwM03w+uvuy3ZdQfXcduU2wjwC2DewHl0i+yW67ScR4AZUTEYfBOffHUWkWIiEubYBzoDW+y1Km+EBoX69jyVgAAoW9atrlqW7l1Kh0kdKBlSkjX3r8mToIAZVmww5AdcEhUR6ev0kH9JRH4QkSYetKs8sFJENgPrgPlKqZ89mJ/H+fPwn8SejrXbjKyJiIAFC+Cee/Kc1Ozts7n9u9upXqo6K+9bSY1SNfKcZsXiFSkSUMSIisHgw7ja/PWyUup7EbkZ6AS8C3wKtPCEUUqpvcCNnkjbLhJTEkm+7JlOcLcREQH79kFcXJ6S+WrTVzww9wGaV2zO/LvnU7pIabeY5yd+1Ayvya6TRlQMBl/F1eYvx7ClbsAEpdR8wPVhOwbfH/0FWlTyOPprzJox3D/nfjpW78iSe5e4TVAc1Cxdk5j4GLemWZA4n3yekb+NpPeM3r79EqMUHDxolrAugLhaUzkoIuOBKOAdEQnGR/tjfBV/P38UiuTLyTkaRutVypfPtagopXj515cZ9fso+tTtw+RekwkOCHazgbpf5aedP5GSmpLGjUthJyU1ha82fcWry1/lyDndL7bh0AZaV25ts2Epuva7fXvabccOOHsWKlSAAwe85sXB4Hlc/Vf2A7oCo5VSp0WkAuCjC677Jv6i/zTnks+5/e3dbURE6DfIlJQcnXY59TKPLniUzzZ+xtDGQ/ms+2eZulvJK5HhkaSkphB7OpYbSt/gkTzyE0op5uycw3+X/pcdJ3bQpnIbPuv2GXdOv5NV/6zyjqhcvAixsbB3r9727Lm6HxOT9kWlQgWoUwfuvVef99VXegnrpk09b6fBK7gqKiHAcgARKQ0kAb96yKYCSa3wWuw5tcf3RQWgbtYuU5w5l3yOIT8N4ftt3/Ns62d5u9PbiIiHDEw7Aqywi8qaA2t4ZvEzrDqwilrhtfix/4/cUesORITrS13P6rjV7s1QKfjnH1i7Vm9//gm7d8OhQ2njFSsGNWrADTfA7bdrEalTB2rXhpIlr8Y7ckSLypIlRlQKEK6Kyp9AZeAUIEBJ4IiIHAUeUEpt9Ix5BYdBjQaxYPcC3x5W7BCVu+92KfrfR/+m7/d9iTkZw7tR7/J066c9aJzGWVRur3m7x/PzRbYf385Lv77ED9t/ICI0gvHdxzOk8ZA0zYFtqrRhYcxClFK5F/lz5/RyCGvXwh9/6E/HkPMiRaBxY+jcGa6/XouIYytbVs95yo6ICKhfH5Yuheeey52NBp/DVVFZDMxUSv0CICKdgT7AV8A4PDQKrCCRb9ZUgWznqiil+Dr6ax5d8CglQkqw5N4ltK/e3gsGQniRcEqFlCp0nfWXUy+zIGYBH6//mEV7FhEaFMpr7V7jyVZPZuidoE3lNnyz+Rt2n9x9ZS2aTDlzBnbu1P0cjm37dti1C1JTdZyaNSEqClq21FuDBhAYmPcL69QJPvsMEhMhJCTv6bmRlNQU1h1cx9K9S6lfrj696vSy26R8gaui0lIp9YDji1JqkYiMVkoNszrtDdnwzqp3AB8XlfLl9edbb8GAARlGOZ98nkcWPMI3m7+hY/WOTLlrCuVDy3vNRBEpVMOKT148yVebvmLc+nHsO72P68KuY2T7kQxrOoxyxcplel6bym0AWHVglRaV1FQ9VNxZOBzb4cNXTwwI0AJSpw707asFpEULCPfQwmgdO8IHH8Dq1dChg2fycBGlFLtP7mbx3sUs3ruYZfuWcTbpLKBfZnrU6mEGh7iAqyV0WESeA6ZZ3/sDR0XEH0j1iGUFjPPJ5wEfF5VSpXSzxalTGQZvPbaVvt/3ZceJHYy4dQQvtX3JYx3yWREZHsmK/fnWY49LbD6ymY/WfcSUv6eQmJJI26pteafTO9xZ+04C/bOpIRw8SJ2Vf1OSEFZPeoPBQ8bqmsiFC1fjlCyphaNrV93X4diqV3dPDcRVbr1Vj/xassQWUUlKSWJ+zHx+3v0zi/YsYv+Z/QBUK1mNAfUGEHV9FAlJCQyZM4Rf9/1K1PVRXrcxv+GqqNwNvAr8aH1fZR3zR48MM2SD8+gvn8XPTz9QMhhWPCl6Eo8seITQoFAW37uYjjU62mCgJrJ0JJP/mszFSxcpEujZtV+8zfLY5bzy6yv8/s/vFAkowr0N7+XR5o/SsHwWTj4vXoTff4dfftHb1q34Aa3+BavK/APlO+qHt7N4lCvnWr+HpwkL07WhpUu9mu3WY1v5ctOXfLP5G+IvxlM8uDgdqnfguTbPEXV9FNeXuv5KX1RiSiKP//w4M7bOMKLiAi6JirVQ1mOZBO92nzkFF8cbfUJSgs2WZENQUBpRSb6czMPzHuar6K9oV60d3931HRXCKtho4NXO+t0nd9OgfANbbXEXm49s5r9L/8vPu3+mUvFKvNf5Pe5rdB+lipS6NrJSsHWrFpBFi2DFCt0nERwMt9wCgwZBhw60OT2XhStf4+QPU3x3xCHofpWRI3UNuVQG1+smziefZ/rW6Xzx5xesiVtDoF8gd9a+k6FNhtKheodMm7ZCAkLoWasnP+z4gXHdxmVfUyzkuCQqIlIWeBaohx5eDIBSyt5G0HxEvqipgBaVpKQrXx9f+DhfRX/FS7e8xKvtXvWJNmXnEWD5XVT2ndrHy7++zHd/f0fJkJK8G/Uu/3fT/6WtgZ0+DevW6e2PP/TmWEOobl146CHo0gXatoWiRa+c1iY2AVa+xpoDeXfm6VE6doTXXoPly6GXezvDlVJsPLyRzzd+ztQtU0lITqB2mdq81/k97m14L2WLlXUpnX71+jHl7yks3beUrjd0dauNBQ1XnxBTgOlAd+AhYBDg+ytj+RCdanRi5YGVvi8q1arpUT/Al39+yfiN43muzXOM7DDSXruccMxPiTmZf0eAHT9/nDdWvMGnGz4lwC+A59o8x3M3P0fJwDA9GXDt2qsisnPn1RPr1IFu3fQSBZ07Q+XKmeZx03U34S/+rD6w2rdFpUULPbdlyRK3icruk7uZ+vdUvtvyHTtO7KBIQBH61+/P0MZDaV25dY6HWXe+vjPFg4szY+sMIyrZ4KqohCulvhSR/yilfgN+E5H1njSsoPFqu1f53+r/+b6o3HYbbN7Mun/W8MiCR4iqEcWoDqPstioNYcFhVAitkC+9FZ9LPseYNWN4d/W7XLx0kfsbD+HVCgO5bvXf8PFg/bZ+5oyOXL68fuDee6/+vOkmKFHC5byKBRWjcYXGrDqwyiPX4jaCgnSfTx77VQ6ePcj0rdOZumUqGw5tQBBuqXoLw1sMZ0D9AZQIcb3s0uNoApu9Yzafdf/Md10t+QCuisol6/OwiHQDDgE+3Ejrm/j8mioAEREcC7lM7xl9uC7sOqb2nmrLCK/syG/r1aeqVCZFT+LFZS9y+Nxheoc0YdSeCGp98hMc+1xHqlFDD+Pt0AFat4YqVfLcmd6mchsmbJzApcuXfLsvoFMnvezCgQNZ1r7Sc+LCCWZtm8XULVNZsX8FCkXTCk0ZHTWa/vX7U6l4JbeZ2K9eP77961uW7F1SaCfeuoKrovKGiJQAngI+AooDwz1lVEHktim3kZCUwLlLvi0ql776gn594cTFeFb/aw3hRT00PyGPRIZH8uOOH+02wyWW71rEk3P+j03nd9PyeDCzfoJWcX9qP1hRUbpPoX173fToZtpUbsOHf3zIpiObaF6xudvTdxsdrdGES5fC4MFZRk1ISuCnnT8xdctUFu1ZREpqCrXL1GZEuxEMrD8w+8meuSSqRhQlgkswY+sMIypZ4Oror3nW7hmgPYCIDPeQTVjpdwU+RA9b/kIp9bYn8/M0Fy9dxE/8fL6m8my9Q/xWDb6t8R8aV2hstzmZEhkeyfELxzmdeJqSISXtNuda4uLYPWciz8R8wo8lj1DlNEz9LZD+FTogz92u38xr1fL4sF6HQ8nVB1b7tqjUr6+HOS9ZkqGoJKYksjBmIVO3TGXurrkkpiRSpUQVnmz5JAMbDOTG8jd61OccQHBAML3q9GL29tkkpSR5xAu3xzlzBrZv5/CWNczb65l1D/MylOdJ4AM32ZEGa1LlJ2hX+3HAehGZo5Ta5on8vIW/n79Pi8p3f3/HBzcc5/G1cE/1enabkyWOEWAx8THcVPEmm61BD/PdsgVmzOD0zz8ysvQWPmoBwcWEUeda8ETH5ygyqkua0VneoGLxilQtUZVVB1YxvOVwr+adI/z8dG1l6VJdliKkqlSW7l3Kd1u+44ftP3A26SzlipXj/sb3M7D+QFpVboWfeHcFjn51+zExeiKL9y6me2R3r+btEhcuwP792mv0vn36MzaW1H172XR+DwvKnmZuLVhfEfCQJuZFVDz5WtAc2G2tAImITAN6AvlbVMTfZ+epRB+JZuicodxyqjijF52FDkftNilLapbWTRy74nfZKyrbtsH06TBjBqk7dzChmfBSZ39OBglDqvXijbs+JsLmeT1tqrTh132/5s25pDfo1AmmTuXM5j/4OnEtn6z/hN0nd1M8uDh31bmLgfUHZjmfxBt0rNGRUiGlmLF1hr2i4ljk7M8/024HD16JElccFkf6s6h+URZ3SSQ+8BKC0KJoJKNq3k6P5vfScIT7V4XPy6+j3GbFtVQEDjh9jyOd00oRGQYMA6hSpYoHTXEfvlpTOXnxJHdNv4tSRUoxY1tVAvkjW6eSdlOjVA38xM+ezvqdO2HGDC0mW7eCCP90acl9/67FsuSdtKt2M+93eZ9GEY28b1sGtKnchu/+/o7Y07FUL1XdbnMyZXvTqnx8O0z66VbOk0ybym0Y2X4kd9a+k5AA33A2GeQfRK/avZi5fSaJKYnesevSJb3EwNatsGnTVQFxrJrp5we1a3Oh/c38ViuYRWHHWHRpB9vOxwKXiQgtRrfre9G5Rmc61ejkcV99WYqKiCSQsXgIYKt/DKXUBGACQLNmzTwpcG6he2R3ziWf48SFE3abkobLqZcZOGsgcWfjWHHfCiJSV8LGnT4vKsEBwVQrWc17jiXj4uC77/S2ebPuC7nlFtRHH/FNg1QeX/MyqSqVCd0nMLTJUJ+qETj3q/iaqDi8L49dN5Yle5cQ3FQYeKIcj73wE00quP8t2h30q9ePr6K/4pfdv9Czdk/3JXz+/FUP0c7b7t1XF84LCIB69fRcpSZNSG3SmOWlTjNp5wxmbZvF+UvnCTkXQtuqbRlS41E6X9+Z+uXqe/V+zFJUlFJh3jIkHQfR67c4qGQdy7c83fppYk/HMm3LtOwje5ERy0ewaM8ixncfT8tKLeHplvDDDz4vKuCFYcVnz+qy+PZb+PVX3eTQsiV8+CH07s2xkoEMmzuMn5b/xC1VbmHinROpUaqG5+zJJQ3KNSAsKIxVB1bxr4b/stscAM4mneWrTV8x9o+x7Du9j4phFRnVYRQPTN1F2W9nwRjf9ZTQoXoHShcpzYxtM3InKidOXCsc27frBdAcBAToRc7q1IG77rq60Fm9ehASQkx8DJM2T+LbdXfzz5l/KB5cnIH1B9K3Xl9uqXKLrT7x7Pe5kTHrgZoiUh0tJgPQDizzNb42T2XJ3iWM+n0U9zW6j2FNh10NiIi4Mqvel4ksHcnKf1a6t6/g0iXtT+vbb+Gnn7RPrRtugFdfhXvu0QtSAT9s/4EHv3uQhKQERkeNZnjL4T45nwd0s2vLSi19YhLkP2f+YewfY/n8z885m3SWNpXbpPW+fGwWfDoJ1q/Xc3V8kED/QO6qfRfTtk7L3KmpUnrOTUbiccKptaJIEe3g8+abrwpH3br6PgtKO8HydOJpZmz9hkmbJ7H6wGr8xI+oGlG83fFt7qx9p884V/VJUVFKpYjIo8Av6CHFXymlttpsVp5oN7Ed+8/sJ+lykk9MRDt67ij3/HAPtcvU5qPbProa0K6dXlfcyf+Xr1IzvCbnks9x9PxRIkIjcp+QUvohNnkyTJum/WqFh8P992shadHiytDf04mneXzh43z717c0qdCEb+78hnrlfHukHOh+ldd+e40ziWfyNLM8t6w/uJ4xa8fw/dbvAehbry9PtHzi2mHO7dvrsl6yxGdFBXQT2BebvuDn3T/Tq8zNej2YbduuCseOHXrlTAelS2vBuPPOq+JRp46e4OqX+Qi22NOxzN81n3kx8/h1368kXU6ibtm6vNPpHe5peA/XhV3n+YvNIT4pKgBKqQXAArvtcCcOp5LnL52npH9J2+xIVancM/seziSdYcm/l1AsqFjaCEFBet3x5ORr3pZ8CWfHkrkSlX37tJBMnqxrZsHB0LOnFpIuXa659mX7ljHox0EcTjjMK21f4aW2L9n+cuAqrSu3RqH44+AfdL6+s1fyvJx6mbm75vLemvdY+c9KigcXZ3jL4TzW/DGqlqya8UmlS0OTJlpUXnnFK3bmmIsXab/nMmVUEWaMGUKvL05fDatUSYvFkCFpxcPFJZZTUlNYG7eWebvmMT9mPluObQH0vf7ITY9wd4O7aVqhqU/12aXHZ0WlIOJoHjmXfM7WCXtvr3ybJXuXMKH7BOqXq39tBMfD9Ngx/SfxUZxFpW3Vtq6ddPKkHrk1eTKsspqD2rXTa6T37p2hb62klCReWvYSo9eMplZ4Ldbcv8Y35sbkgJaVWuInfqz6Z5XHReVM4hkmRk/ko3UfsefUHqqWqMqYzmO4v8n9FA8unn0CnTrBe+/pN/3Qa5dK9jqXL2snn4sXa7FbuZKApCR63+HH5BuTuDDyFYp27KoncIblvBs6MSWR+bvmM3vHbBbuXsjJiycJ8AugbdW2jOk8hm6R3a7c6/kBIypexBfc36/6ZxWv/PoK/ev1Z2iToRlHcqz8d+SIT4tK5eKVCfYPzr6zPiVFrz3y5Zcwb57uN6lbVy+bfPfdugkiE7Yf387dP9xN9JFoHmr6EO91eY+igd6dwOgOwoLDaFi+oUf7VXbF7+KjPz5i4uaJnEs+R+vKrXmz45vcVeeunM0t6dQJ3nlHLzx2220es/caHHM/tmxJu23bphdCA2jYEP7v/yAqin5VLzN+RncW3tWQ3nVb5SirVJXKiv0rmPzXZGZum8mZpDOUKVqGHpE96B7ZXbuEsaGZ0h0YUfEidi/UFX8hnoGzBlKtZDUm9JiQeRXaUVM56tsTIP39/Lmh9A2Zi8q+ffDVV/D11/phUa4cPPaYbt5q1CjL5gilFJ9u+JSnFj1FaFAocwbMoUetHp65EC/RpnIbJkZPJCU1xW0TCFNVKov2LGLsH2NZuHshgX6BDKg/gMdbPE6z65rl0tA2uilyyRLPisqhQ/Dbb1q8/vpLC4jDQzRo32z16+v1apo10zP+y1+d49E2NYVyxcoxY9sMetft7VKWfx39i8l/TWbqlqnEnY0jNCiUu+rcxT0N7qF99fY+sV5RXsn/V5BP6FevH7vid7Hl2BZbaipKKYbMGcKRc0dYff/qzJsh+vXTTUSbNuWbYcU7453WG0lKgh9/hC++0A8lPz+9DvtHH0H37i6tv37s/DGG/DSE+THz6XpDV77u+XXeBgL4CG0qt+GT9Z/w99G/8+zXLTElkS///JKP1n3EzvidRIRG8Fq71xjWdFjey6pIES0sS5bkLZ30HDyolxb47Tf9GWOtx1O8uH7J+Ne/tIjUr6+H7pbO2hF7gF8Avev0ZtLmSZxPPn9t36TFnpN7+H7b90z5ewpbjm0hwC+Arjd05d2od7mj1h35suabJUqpfL81bdpU5Qc2HNygGIGas2OO1/P+YM0HihGo99e8n33kixeVAqXeeMPjduWVZxc9q4JGBqmUg3FKDR+uVHi4tr1qVaVef12pf/7JUXrzd81X5d4tp4JHBquxa8eq1NRUzxhuA7GnYhUjUB/98VGe0pm3c566/sPrFSNQN024SU3ePFklpSS5yUqLN9/Uv+PRo7lP4+hRpb77TqmhQ5W64QadHihVooRSPXoo9d57Sm3YoFRKSq6z+HXfr4oRqOlbpqc5HhMfo95c8aZq/FljxQgUI1CtvmilPln3iTp+/njur8nNABuUm5/HpqbiJS5cumBbn8qGQxt4ZvEz9IjswX9a/CfryBcu6M8SJfJNTSX5cjL/tKlH9YMX9JDNoUN1u3wWQzXTk3w5mWcXP8uHf3xIw/INWfrvpRkPYsjHVClRhYphFVl1YBWPNn80x+fvPbWX4T8PZ+6uudQuU5tF9ywi6vooD1iK/v1eeAGWLYMBA1w758IFWLlSd6gvXqw9HwCULKmXWn7kET0oo2FD8HfPnKJbqtxC+WLlmbF1Bo0iGvH91u/5ftv3bD6q825RsQWjo0bTp26fzEe8FTCMqHiJ26fcTtJlPffDm6JyNuksA2YOoHxoeb7u+XX2QxFvt9aJiIjw+T4V4uOJ/Pg7qAq76pan+sKf9ESyHBJ3No6+3/dlbdxaHm/+OO9EveMzvqbciYjQpkobVv2Ts876i5cu8s6qd3h75dsE+AXwv07/4z8t/+PZ1Q+bNNFisGRJ5qKSmKhHZS1frkVk1Srd/BkYqJvPRo3S69U0aeI2EUmPv58/fer24ZP1nzBr+ywAWlVqxZjOY+hdtzdVSuQPv4TuxIiKF7GjpvLQvIeIPR3L8sHLc7bgVkSEb9dU5s+HoUOJTDwBw2HXSw/TJReCsnjPYu7+4W4SUxKZ0WcGfev1db+tPkSbym2YsXUGB84coHKJrFdYVEoxd9dchv88nH2n9zGg/gBGR42mYvGKnjfU319PhFy8WDdaXboEf/8NGzZc3bZsueoTq0GDK6OyuOUWvea9l3is+WPEnY2jXbV29K7TO9tyLegYUfEizvNUvMGCmAVM3TKV19q9xs1Vbs7ZyRER2hOqr3H2LDzxhB7V1bAh5SYtpPjPt7Lr1J4cJZOqUhm1YhSvLn+VumXrMqvfLGqVqeUho30HZ+eS/Uv0zzTervhdPPHLEyyIWUDdsnVZ9u9ltK/e3ltmajp1gtmzdU1j2zY9GRd0B3qzZvDss/qzVSt9v9pErTK1+HHAj7bl72sYUfEighASEOIVUUlMSeTxhY9TK7wW/735vzlPwBdrKsuWwX33aY/BL7wAr7yCBAdTc11NYk7GuJxM/IV47p19Lwt3L+SehvfwWbfPMh25U9C4sfyNFA0syqoDq+hf/1pR2XR4E2+tfIuZ22YSGhTKmM5jeLT5o/Z4DrjjDj1qLzwchg/XAtKsmV522YdnlBd2jKh4mdCgUBKSPT9PZfTq0ew5tYdF9yzKXdt3+fKQkKA7P728WuE1bNsGo0fr+SaRkdrPUoury+tEhkeyNm6tS0mtP7iePt/34ci5I3za7VMebPqgT7u8cDeB/oG0qNgizSRIpRQr9q/grZVv8cueXygeXJzn2jzH8JbDPb72RpZUqqT9aBnyFUZUvMTgRoMBeO231zxeU9l/ej9v/v4mvev0zvnoHMf64MpaouboUahuwxocSulO2jFj4OefISQEnnwSRo68RuQiwyOZtmValuuGK6UYv3E8//n5P1QIrcCqIatyPzkvn9OmchveWvkWCUkJLI9dzlsr32JN3BrKFSvHWx3f4uFmD+fb2dwG+zGi4iUcojJmzRiPi8oTvzyBiDCmy5icn+wQlYUL9eeRI94VlcREmDIFPvhAd8SWL6+F5MEHtVO+DIgMj0Sh2HNqD3XL1r0m/OKlizw8/2EmbZ7EbTfcxre9vs3ZoIUCRuvKrbmsLlPnkzocTDhItZLV+OT2T7iv0X0+4z7dkH8xouIlHCs+enpNlV92/8LsHbMZ1WFU7oYzOtZ6cHR8eqtf5ehR+PRTGDdOu55v2BAmTtTDSYMzrn04cHYsmV5UYk/Hctf0u9h0ZBMjbh3By7e+jJ+4Pn+lINK6cmtKBJegZEhJ3u70NgPqDygQ7kEMvoG5k7xEnxl9AM+KSlJKEo8tfIyapWvyVKuncpdIH20nU6fqT0+LilJ6JNdjj2mnfd266WYux7oaLlCzdE0AYuLTdtYv2rOIgbMGahfsA+fSPbK7283Pj5QIKcGhpw4REhBS6AXW4H587o4SkREiclBEoq3tdrttcieeFJUxa8YQczKGsbeNzbRvwWUc6z94UlQSEuDee/UM+NatdafsvHnQoUOORveUCClBuWLlrjiWVErx1u9v0XVyV64Lu44NwzYYQUlH0cCiRlAMHsFXayrvK6VG222EJ/CUqPxz5h/e+P0N7qx9J11v6Jr3BAMCoEwZz82q37xZO6/cvRtef10PEc7DrOfI8Eh2ndzF2aSzDP5xMLN3zGZA/QF80eOLQjNc2GDwBXxVVAosnhKVpxY9RapK5f0u77svUU/MVVEKJkyA//xHT2JbtgxuvTXPyUaWjmT2jtk0/7w5u0/uZkznMQxvObxQDRc2GHwBX63/Pioif4nIVyJSKqMIIjJMRDaIyIbjx497275c4wlRWbJ3CTO3zeSFm1+gWslq7kvY3aJy9iwMHKjXp7j1Vu23yQ2CArqmcirxFKcST7Hk30t4otUTRlAMBhuwpaYiIkuAjPwqvAh8CowElPX5HjAkfUSl1ARgAkCzZs2Ux4x1Ew83exiAnfE7uZhy0W0LJSVfTubRBY9yfanreabNM3lOj4cfvrofEaHXbncHmzbp5q59++DNN/XyvTnwIpwdA+oPYP+Z/bxwywtUKu67q1UaDAUdW0RFKdXJlXgi8jkwz8PmeAWHS4z31+jmqfPJ590yweyDtR+wM34n8wbOc49n3f5OrjvKl9d9KkrlzS3GpEkwbJju/F++HG7OoR8yF6hasirjuo1ze7oGgyFn+Fzzl4hUcPraC9hily3u5MCZAxw4c4DQoFDAPU4lD549yOu/vU6PyB50i+yW5/QAOHBAb6BrKomJutkqNyil1xofPFh7jo2O9oigGAwG38EXO+r/JyKN0M1fscCDtlrjJu6dfS8ADzbVl+MOUflg7QckpiTyQdcP8pzWFe7VdrJ8edoJkCVyWKtKTYWnn4b339cTGCdNgiAPrr9hMBh8Ap8TFaXUvXbb4EncVVM5n3yeLzZ9Qe+6valRqoY7TLsWZ1GplQO38JcuwZAhMHmyntT4wQdu7T8xGAy+i8+JSkHHXaIy5e8pnE48zePNH3eHWRlT3vJQm5MRYOfP61n5P/+sV957/nnjptxgKEQYUfEy7hAVpRRj/xhLkwpNriy65BEcNRVXJ0DGx2s3K+vXw+ef65nyBoOhUGFExcu4Q1R+jf2Vrce3MrHnRM/OxShdWs+sd6Wm8s8/0KWLHjI8axbceafn7DIYDD6LERUv4XDw6BCVvCzUNfaPsZQtWjbDlfvyzFNOjij9/HQTWHaism2bFpSEBFi0CNq2db9dBoMhX2BExUv0qNUDgNOJp4Hc11T2ndrHnJ1zeOGWF9wzLyU9PXqk/Z7VrHql4Ntv4dFHoVgxWLFCu6w3GAyFFjMkx0vsPLGTnSd2UixQOzfMraiMWz8OP/G7MkPf7ezcqTcHmdVUTp7UEyUHDYLGjeGPP4ygGAwGU1PxFg/O0/NTlg9eTrB/cK5ExTGMuE/dPlQsXtHdJmoetKYFLV+uPyMi9KRFZ5Yu1WJy7Bi8/baej5IHD8MGg6HgYGoqNpBbp5KT/5qshxG38OAw4vREROjRX6mpkJSkBaRTJwgLg7VrtQ8vIygGg8HCiIoN5EZUlFKMXTeWphWa0qpSKw9ZlgEREXD5Mvz2G9x0E7z3HjzyCGzcCE2aeM8Og8GQLzDNXzaQG1FZtm8Z245vY9Kdk7zr0t0xV6VjR+0Qcv58uL1ALcZpMBjciKmp2EBuRGXsOmsYcT0PDCPOishI/dm9O/z9txEUg8GQJaam4iVeavvSlf2cisq+U/uYu3MuL97yYt7Xns+Ol15K+/3GGyEuDq67zrhbMRgM2WJExUt0qnF1CZmw4DCOX3B9tcpP1n+Cv58/DzV7yBOmpaVTBkvdVPTQSDODwVDgMM1fXiL6SDTRR6KBnNVUziWf44s/PTyM2Jno6GuHEBsMBoOLmJqKlxj+83BAz1MJDXRdVCb/NZkzSWc8643YmeHD9adjnorBYDDkAFtqKiLSV0S2ikiqiDRLF/a8iOwWkZ0i0sUO+zyNqzUVpRQfrfuIphWa0rJSSy9YZjAYDHnDruavLcBdwArngyJSFxgA1AO6AuNEpMDNrAsNCuXCpQtcTr2cZTzHMOLHWzzu3WHEBoPBkEtsERWl1Hal1M4MgnoC05RSSUqpfcBuoLl3rfM8Dk/FFy5dyDLeuA3j7BlGbDAYDLnE1zrqKwIHnL7HWceuQUSGicgGEdlw/LjrI6l8AVfWVElMSeTn3T/Tt25fzw8jNhgMBjfhsY56EVkCRGQQ9KJS6qe8pq+UmgBMAGjWrJnKa3qe5s2Ob17Zd0VUfov9jQuXLtAtspvHbUvDm29mH8dgMBgywWOiopTKYMJDthwEKjt9r2Qdy/c4L/vrykJd82PmUySgCO2rtfe4bWlo7cHliQ0GQ4HH15q/5gADRCRYRKoDNYF1NtvkFlYfWM3qA6sBPfkRMq+pKKWYt2seHWt0pEhgEa/ZCMDq1XozGAyGXGDLPBUR6QV8BJQF5otItFKqi1Jqq4jMALYBKcD/KaWyHiKVT3hh6QuANU8lm+avHSd2sO/0Pp5t86zX7LvCC9pOM0/FYDDkBltERSk1G5idSdgoYJR3LfIu2YnK/Jj5ANxe0zhvNBgM+Qtfa/4qFLgiKg3KNaBKiSreNMtgMBjyjBEVG8hKVE4nnmblPyvpHtnd22YZDAZDnjGiYgNZicqiPYtISU2hW00vDyU2GAwGN2AcSnqJD7p+cGU/yD+IQL/ADEVlfsx8ShcpbZ+vrw8+sCdfg8FQIDCi4iUaRTRK8z00KJSEpLTzVFJVKgtjFtL1hq74+9nk8qxRI3vyNRgMBQIjKl5iyd4lwNXFusKCwzh3KW1NZf3B9Ry/cNzepq8l2s4MF+syGNzEpUuXiIuLIzEx0W5TCgUhISFUqlSJwMBAj+dlRMVLvLHiDeCqqGTk/n7ernn4iR9db+jqdfuu8Ia204iKwZPExcURFhZGtWrVjAduD6OUIj4+nri4OKpXr+7x/ExHvU1kJCrzY+bTunJrShcpbZNVBoN3SExMJDw83AiKFxARwsPDvVYrNKJiE+lF5VDCITYd2WRGfRkKDUZQvIc3y9qIik2kF5UFMQsAjKgYDIZ8jREVm0gvKvN2zaNy8crUL1ffRqsMhsJBYmIizZs358Ybb6RevXq8+uqrGcabOHEiZcuWpVGjRjRq1Igvvvgiy3RnzZqFiLBhw4YMw6tVq0aDBg1o1KgRzZo1yzDO4MGDmTlzpkvX0b59e3755Zc0xz744AMefvhhl873BKaj3kuM7z4+zffQwKuikpSSxJK9S/j3jf+2v0lg/Pjs4xgM+Zzg4GCWLVtGaGgoly5d4uabb+a2226jZctr54f179+fjz/+ONs0ExIS+PDDD2nRokWW8X799VfKlCmTa9udGThwINOmTaNLly5Xjk2bNo3//e9/bkk/NxhR8RK1ytRK8925pvLb/t84f+m8bzR91aqVfRyDwZ0MHw7R0e5Ns1GjLCfyigihodqzxaVLl7h06VKeX+hefvllnnvuOd599908pZM+zQMHDvDll18yZswYZsyYQVJSEr169eK1116jT58+vPTSSyQnJxMUFERsbCyHDh3illtucZsNOcU0f3mJuTvnMnfn3Cvfw4LDOJd8jlSVyvxd8wkJCKF9dS8vyJURc+fqzWAo4Fy+fJlGjRpRrlw5oqKiMq1hzJo1i4YNG9KnTx8OHDiQYZw///yTAwcO0K1b1i+GIkLnzp1p2rQpEyZMyDLuM888w/Hjx/n6669ZunQpMTExrFu3jujoaDZu3MiKFSsoXbo0zZs3Z+HChYCupfTr18/WFg9TU/ES7615D4AetXoAV/1/nU8+z7yYeXSo3oGigUVts+8K72k76dHDXjsMhQebXAP5+/sTHR3N6dOn6dWrF1u2bKF+/bR9mj169GDgwIEEBwczfvx4Bg0axLJly9LESU1N5cknn2TixInZ5rly5UoqVqzIsWPHiIqKonbt2rRt2/aaeCNHjqRFixZXhGfRokUsWrSIxo0bA3Du3DliYmJo27btlSawnj17Mm3aNL788stcloh7sKWmIiJ9RWSriKSKSDOn49VE5KKIRFvbZ3bY5w0corLx8Eb2ntrrG01fBkMhpGTJkrRv356ff/75mrDw8HCCg4MBGDp0KBs3bgTgxRdfvNJ5n5CQwJYtW2jXrh3VqlVj7dq13HHHHRl21lesWBGAcuXK0atXL9aty3hh25tuuomNGzdy8uRJQE9gfP7554mOjiY6Oprdu3dz//33A9CzZ0+WLl3Kn3/+yYULF2jatGneCyUP2NX8tQW4C1iRQdgepVQja3vIy3Z5DYeoTN8yHTBDiQ0Gb3L8+HFOnz4NwMWLF1m8eDG1a9e+Jt7hw4ev7M+ZM4c6deoAMGrUqCsP+BIlSnDixAliY2OJjY2lZcuWzJkz55rRXefPnychIeHK/qJFi66pGTno2rUr//3vf+nWrRsJCQl06dKFr776inPndD/swYMHOXbsGAChoaG0b9+eIUOGMHDgwLwVjBuwa+XH7VC4Jz85RGXm9pnUL1efqiWr2myRwVB4OHz4MIMGDeLy5cukpqbSr18/unfXaxi98sorNGvWjDvuuIOxY8cyZ84cAgICKF26tEtNXM4cOnSIoUOHsmDBAo4ePUqvXr0ASElJ4e6776Zr18xdMvXt25eEhATuuOMOFixYwN13302rVq0ALSSTJ0+mXLlygB4F1qtXL6ZNm5aL0nAvopSyL3OR5cDTSqkN1vdqwFZgF3AWeEkp9Xt26TRr1kxlNi7cV2g3sR2g16gH7WAy6tsoAJ5r8xxvd3rbJsvS0a6d/jRr1Bs8yPbt26+89Ru8Q0ZlLiIblVIZT5jJJR6rqYjIEiAig6AXlVI/ZXLaYaCKUipeRJoCP4pIPaXU2QzSHwYMA6hSxfeX3f2217dpvjtqKuBjTV/ffpt9HIPBYMgEj4mKUirHbm6VUklAkrW/UUT2AJHANdUQpdQEYALomkrerPU8lUtUTvPdISqlQkrRqnIrO0zKmMqVs49jMBgMmeBT81REpKyI+Fv7NYCawF57rXIP07dMv9IpDxAWFAZA1xu6EuDnQyO7p0/Xm8FgMOQCW55mItIL+AgoC8wXkWilVBegLfC6iFwCUoGHlFIn7bDR3Xy64VMA+tfvD0D50PI0qdCE+xvfb6dZ1/KptpP+/e21w2Aw5EvsGv01G5idwfFZwCzvW+R9QgJC2Dhso91mGAwGg1vxqeYvg8FgMORvjKgYDIZCy+XLl2ncuPGVOSrpcdX1vavx2rVrR61ata7Ec0xgdGbEiBGMHj3aJfvvu+8+xqfzLP7jjz9y2223uXS+J/ChHmKDwWDwLh9++CF16tTh7NlrZi1cwVXX967GmzJlSqZrqeSUgQMH8tZbb/Hggw9eOTZt2jRbZ9YbUfESM/u5tuiO7bi4OJDB4C6G/zyc6CPRbk2zUUQjPuj6QZZx4uLimD9/Pi+++CJjxoxxa/7u4vPPP+eHH37ghx9+YNasWYwdO5bk5GRatGjBuHHj6NixI4MGDeLw4cNUqFCB8+fPs2TJkmw9IHsS0/zlJcoULUOZou5ZmMejlCmjN4OhgDN8+HD+97//4eeX9WPQFdf3OYl333330ahRI0aOHElWHk0+/vhj5s2bx48//khsbCzTp09n1apVREdH4+/vz5QpU/D396d3797MmDEDgLlz59KuXTuKFy+ezdV7EKVUvt+aNm2qfJ2vN32tvt70td1mZM/XX+vNYPAg27ZtszX/uXPnqocfflgppdSvv/6qunXrlmG8EydOqMTERKWUUp999plq3759nuLFxcUppZQ6e/asioqKUpMmTbomzquvvqoaNGigbr/9dpWcnKyUUuqjjz5SFSpUUDfeeKO68cYbVWRkpHr11VeVUkqtXLlStWzZUimlVM+ePdXMmTMzzDujMgc2KDc/j01NxUtMjJ7IxOiJdpuRPRMn6s1gKMCsWrWKOXPmUK1aNQYMGMCyZcu45557ronniuv7rOKlx+H6PiwsjLvvvjtT1/cNGjQgNjaWuLg4QL/8Dxo06Ipn5J07dzJixAgAWrduzeHDh9m8eTOrV6/OdqEwT2NExWAwFDreeust4uLiiI2NZdq0aXTo0IHJkydfE88V1/dZxXMmJSWFEydOAHoJ43nz5mXq+r5x48aMHz+eO+64g0OHDtGxY0dmzpx5ZbTYyZMn2b9/P6C9vffv359BgwZx2223ERISkosScR9GVAwGg8GJV155hTlz5gAwduxY6tWrx4033sjYsWMzdX2fVTxHbSYpKYkuXbrQsGFDGjVqRMWKFXnggQcytePmm29m9OjRdOvWjXLlyvHGG2/QuXNnGjZsSFRUVBohGzhwIJs3b/aJ9VRsdX3vLvKj63ufxbi+N3gB4/re+3jL9b2pqRgMBoPBbZh5Kl5iwb8W2G2CayzIJ3YaDAafxIiKlygaWNRuE1yjaD6x05DvUUoV6iXFvYk3uzlM85eXGLd+HOPWj7PbjOwZN05vBoMHCQkJIT4+3qsPu8KKUor4+HivjQozNRUvMWOrnvH6yE2P2GxJNlgzc3nEx+005GsqVapEXFwcx48ft9uUQkFISAiVKlXySl5GVAwGg9cJDAykevXqdpth8AC2NH+JyLsiskNE/hKR2SJS0inseRHZLSI7RaSLHfYZDAaDIXfY1aeyGKivlGoI7AKeBxCRusAAoB7QFRjnWLPeYDAYDL6PLaKilFqklEqxvq4FHI19PYFpSqkkpdQ+YDfQ3A4bDQaDwZBzfKFPZQgw3dqviBYZB3HWsWsQkWHAMOtrkohs8ZiF7qOM3Ccn7DbCBcog+cROMHa6j/xgZ36wEfKPnbXcnaDHREVElgARGQS9qJT6yYrzIpACTMlp+kqpCcAEK50N7nY14AmMne7F2Ole8oOd+cFGyF92ujtNj4mKUqpTVuEiMhjoDnRUVwerHwQqO0WrZB0zGAwGQz7ArtFfXYFngTuUUhecguYAA0QkWESqAzWBjBccMBgMBoPPYVefysdAMLDYctOwVin1kFJqq4jMALahm8X+Tyl12YX07FuQOWcYO92LsdO95Ac784ONUIjtLBCu7w0Gg8HgGxjfXwaDwWBwG0ZUDAaDweA2fFJURKSr5aZlt4j8N4PwYBGZboX/ISLVnMIydPOSXZretFNEokRko4j8bX12cDpnuZVmtLWVs8nGaiJy0cmOz5zOaWrZvltExoob/Jfnwc5/OdkYLSKpItLICnNrWbpoZ1sR+VNEUkSkT7qwQSISY22DnI7bUZ4Z2ikijURkjYhsFe1Gqb9T2EQR2edUno3sstMKu+xkyxyn49Wte2S3dc8E2WWniLRPd38misidVpgd5fmkiGyzftulIlLVKcw996dSyqc2wB/YA9QAgoDNQN10cR4BPrP2BwDTrf26VvxgoLqVjr8raXrZzsbAddZ+feCg0znLgWY+UJbVgC2ZpLsOaAkIsBC4zS4708VpAOzxRFnmwM5qQEPgG6CP0/HSwF7rs5S1X8rG8szMzkigprV/HXAYKGl9n+gc187ytMLOZZLuDGCAtf8Z8LCddqa7B04CRW0sz/ZO+T/M1f+72+5PX6ypNAd2K6X2KqWSgWlo9y3O9AQmWfszgY6Wembm5sWVNL1mp1Jqk1LqkHV8K1BERILzaI9bbcwsQRGpABRXSq1V+o77BrjTR+wcaJ3rKbK1UykVq5T6C0hNd24XYLFS6qRS6hTa/11Xu8ozMzuVUruUUjHW/iHgGFA2j/a43c7MsO6JDuh7BPQ9c6eP2NkHWKjSTqNwJ67Y+atT/s4ustx2f/qiqFQEDjh9z8hVy5U4SvsQOwOEZ3GuK2l6005negN/KqWSnI59bVWHX85jU0hebawuIptE5DcRucUpflw2aXrbTgf9ganpjrmrLF21M6fn2lWe2SIizdFvvHucDo+ymk7ed8OLUF7tDBGRDSKy1tGkhL4nTqurvgV9pjzRNez096ed5Xk/uuaR1bk5vj99UVQKDSJSD3gHeNDp8L+UUg2AW6ztXjtsQzd7VFFKNQaeBL4TkeI22ZItItICuKCUcvYB5ytlme+w3lC/Be5TSjnevp8HagM3oZtJnrPJPAdVlXaFcjfwgYhcb7M9mWKVZwPgF6fDtpWniNwDNAPedXfavigqrrhquRJHRAKAEkB8Fud6wv1LXuxERCoBs4F/K6WuvAkqpQ5anwnAd+TNS3OubbSaEOMtWzai31YjrfjOS8jZXpYW17wFurksXbUzp+faVZ6ZYr08zEf76bvi4FUpdVhpkoCvsbc8nX/fvej+s8boe6KkdY/kOE1P2GnRD5itlLrkOGBXeYpIJ+BFtEeTpGzOzfn96a5OIndt6Fn+e9Ed7Y7Opnrp4vwfaTttZ1j79UjbUb8X3XmVbZpetrOkFf+uDNIsY+0HotuFH7LJxrKAv7Vfw7qRSquMO+5ut6ssre9+ln01PFWWrtrpFHci13bU70N3gpay9m0rzyzsDAKWAsMziFvB+hTgA+BtG+0sBQRb+2WAGKxOaeB70nbUP2KXnU7H1wLt7S5PtPDuwRqM4Yn7M9cX4MkNuB29eNce9NsSwOtoZQUIsW6c3dYFOz9MXrTO24nTKIWM0rTLTuAl4DwQ7bSVA4oBG4G/0B34H2I92G2wsbdlQzTwJ9DDKc1mwBYrzY+xPDPY+Ju3Q7v6cU7P7WXpop03odudz6Pfmrc6nTvEsn83ulnJzvLM0E7gHuBSunuzkRW2DPjbsnUyEGqjna0tWzZbn/c7pVnDukd2W/dMsM2/ezX0S49fujTtKM8lwFGn33aOu+9P46bFYDAYDG7DF/tUDAaDwZBPMaJiMBgMBrdhRMVgMBgMbsOIisFgMBjchhEVg8FgMLgNIyqGfE06T7XR4uSx2lcRkWYiMtYN6YiILHN4OhCRc+nCB4vIx1mc311EXs+rHQaDM3YtJ2wwuIuLSqlGGQVYvr5EXXUz4hMopTYAG9yQ1O3AZqXU2VyePx8YKSJvK885OTQUMkxNxVCgEL0OzE4R+QY9YauyiDwjIustx32vOcV9UUR2ichKEZkqIk9bx5eLSDNrv4yIxFr7/iLyrlNaD1rH21nnzBSRHSIyxeG8UkRuEpHVIrJZRNaJSJgVf54VXkxEvrLCNolIT+t4PetYtJVXzQwu91/ATy6Wi3Nt7qKI3Kr0JLXlQPfclLXBkBGmpmLI7xQRkWhrfx/wBFATGKSUWisina3vzdFuJuaISFv0zOcBQCP0/+BP9Az8rLgfOKOUusnyKLtKRBZZYY3RboIOAauANiKyDpgO9FdKrbeaqS6mS/NFYJlSaoiIlATWicgS4CHgQ6XUFNGLTPlnYE8b0jojdS4L0K435gA4anMi0gN4FlhtxdmAdrY5I5trNxhcwoiKIb+TpvnL6lPZr646QuxsbZus76FokQlDO/i7YJ13ZeXALOgMNJSrK/uVsNJKBtYppeKstKLRrjnOAIeVUusBHM1U6TzwdwbucNSS0O5oqgBrgBctx6M/KGuNk3SUVtpZZmZlMRjtYsPxvSbaK217ddWx4TH0YlwGg1swomIoiJx32hfgLaXUeOcIIjI8i/NTuNo0HJIurceUUs7uyxGRdoDzejiXcf2/JUBvpdTOdMe3i8gfQDdggYg8qJRalt5OEfFzpc9IRELRtZEHlFKHnYJCuLb2ZDDkGtOnYijo/AIMsR6qiEhF0WvVrwDuFJEiIhIG9HA6JxZoau33SZfWwyISaKUVKSLFssh7J1BBRG6y4oc5uWR3TvMxpz6YxtZnDWCvUmosut+kYSbp18jy6q/yFfC1Uur3dMcj0X1PBoNbMKJiKNAopRah11JZIyJ/o13ghyml/kT3d2xGu/Ne73TaaLR4bEK7VXfwBbAN+FNEtgDjyaJGovSSrv2Bj0RkM3qJ1pB00UaiXfP/JSJbre+g19/YYjWl1Ucv45qe+WgPzVkiIlXR4jjEqbPe0SzW3krHYHALxkuxwQCIyAjgnFJqtN22uIro1QS/UUpF5fL88sB3SqmO7rXMUJgxNRWDIZ9i9Y18Lrlf5rkK8JQbTTIYTE3FYDAYDO7D1FQMBoPB4DaMqBgMBoPBbRhRMRgMBoPbMKJiMBgMBrdhRMVgMBgMbuP/AV0TwObWCJV1AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.figure()\n", "plots = []\n", "colors = ['r','g']\n", "energies = ['3.5-4.5 keV', '4.5-5.5 keV']\n", "\n", "# Plot lag-frequency spectrum\n", "for i in range(0,len(lags)):\n", " plots += plt.plot(cross[i].freq, lags[i], colors[i], label=energies[i])\n", " plt.axvline(v_cuts[i],color=colors[i],linestyle='--')\n", " plt.axhline(h_cuts[i], color=colors[i], linestyle='-.')\n", "\n", "# Define axes and add labels\n", "plt.axis([0,0.2,-20,20])\n", "plt.legend()\n", "plt.xlabel('Frequencies (Hz)')\n", "plt.ylabel('Lags')\n", "plt.ylim(None, 25)\n", "plt.title('Energy Dependent Frequency-lag Spectrum')\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.7" } }, "nbformat": 4, "nbformat_minor": 1 }